JP2003028819A - Method for calculating heat transfer by piston effect at supercritical fluid - Google Patents

Abstract

PROBLEM TO BE SOLVED: To provide a method for calculating a heat transfer by a piston effect at a supercritical fluid in which the heat transfer can be highly accurately and easily calculated with the use of an information processing device such as a computer while a change of physical properties values by a temperature and a pressure is taken into account. SOLUTION: The method for calculating the heat transfer of the supercritical fluid consists of a first step of inputting an analysis region, an initial condition and a criticality, a second step of obtaining physical properties values on the basis of an initial temperature and an initial pressure, a third step of calculating a temperature T with the use of physical properties values obtained in the second step, a fourth step of calculating a pressure p with the use of physical properties values obtained in the second or a sixth step, a fifth step of calculating a heating value Q with the use of physical properties values obtained in the second or the sixth step, and a seventh step of calculating the temperature T with the piston effect taken into account with the use of the heating value Q and physical properties values obtained in the sixth step.

Description

【発明の詳細な説明】Detailed Description of the Invention

【０００１】[0001]

【発明の属する技術分野】本発明は、超臨界流体におけ
るピストン効果による熱移動の算出方法に関する。TECHNICAL FIELD The present invention relates to a method for calculating heat transfer due to a piston effect in a supercritical fluid.

【０００２】[0002]

【従来の技術】近年、超臨界流体の多方面での利用が期
待されている。超臨界流体は、粘度、拡散速度等の性質
が気体と液体の中間に当たる状態にあるため、気体のよ
うにどこにでも拡散する性質と液体のように成分を溶解
する性質とを併せ持ち、例えば、食品や医療品等の分野
では、有機媒体の代替物として利用することができ、環
境に優しい技術として注目されている。ところで、超臨
界流体を工業的に利用するためには、超臨界流体におけ
る熱に関する、例えば、超臨界流体内での熱移動を正確
に把握することが不可欠とされている。2. Description of the Related Art In recent years, the use of supercritical fluids in various fields is expected. Since supercritical fluids have properties such as viscosity and diffusion speed that are in the middle of gas and liquid, they have the property of diffusing anywhere like gas and the property of dissolving components like liquid. It can be used as a substitute for an organic medium in the fields of medicine, medical products, and the like, and is attracting attention as an environment-friendly technology. By the way, in order to industrially use the supercritical fluid, it is indispensable to accurately grasp the heat transfer in the supercritical fluid, for example, the heat transfer in the supercritical fluid.

【０００３】従来、超臨界流体中では、熱拡散率は微小
であるため、自然対流が発生しない場合には、熱移動は
ほとんど無いものと考えられてきた。ところが、ドイツ
のＮｉｔｓｃｈｅとＳｔｒａｕｂ等の１９８６年の報告
によれば、宇宙実験において、対流の発生しない微小重
力下において超臨界流体中で高速な熱移動現象が見られ
ることが明らかになった。この超臨界流体における熱移
動現象は、フランスのＺａｐｐｏｌｉ等によって、１９
９０年にピストン効果と名付けられた。Conventionally, since the thermal diffusivity is very small in a supercritical fluid, it has been considered that there is almost no heat transfer when natural convection does not occur. However, according to a report by Nitsche and Straub et al. In Germany in 1986, it was revealed in space experiments that a high-speed heat transfer phenomenon was observed in a supercritical fluid under microgravity in which convection did not occur. The heat transfer phenomenon in this supercritical fluid was reported by Zappoli et al.
It was named the piston effect in 1990.

【０００４】ピストン効果とは、上述のごとく、超臨界
流体の臨界点近傍における高速な熱移動現象を指し、従
来の熱移動現象である、熱伝導、対流、輻射とは異なる
伝熱機構とも言われているが、理論的には未だ不明な点
も多い。このピストン効果は、前出のＮｉｔｓｃｈｅ、
Ｓｔｒａｕｂ、および日本の小貫等によってさらに研究
が進められ、特に、小貫等はピストン効果に関係ある物
性値を一定とした温度分布に関する基礎式を基に熱移動
現象を解析している。As described above, the piston effect refers to a high-speed heat transfer phenomenon in the vicinity of the critical point of a supercritical fluid, and is also called a heat transfer mechanism different from the conventional heat transfer phenomena of heat conduction, convection, and radiation. However, there are many unclear points in theory. This piston effect is due to Nitsche,
Further research has been conducted by Straub and Konuki of Japan. In particular, Konuki et al. Analyze the heat transfer phenomenon based on a basic equation relating to temperature distribution with constant physical properties related to the piston effect.

【０００５】[0005]

【発明が解決しようとする課題】しかし、上記従来の超
臨界流体のピストン効果による熱移動の解析は、未だ研
究段階にとどまるものであり、工業化レベルでの利用が
有効な算出方法を提供するものではない。また、上記従
来の研究では、ピストン効果による熱移動の算出に必要
な物性値を一定値として計算しているために、実際の超
臨界流体の熱移動現象に即した精度良い算出結果が得ら
れないという問題もあった。However, the conventional analysis of heat transfer due to the piston effect of the supercritical fluid is still in the research stage, and it provides a calculation method that is effectively used at the industrial level. is not. Further, in the above-mentioned conventional research, since the physical property value necessary for calculating the heat transfer due to the piston effect is calculated as a constant value, an accurate calculation result according to the actual heat transfer phenomenon of the supercritical fluid can be obtained. There was also the problem of not having it.

【０００６】本発明の目的は、上記の従来の問題点に鑑
み、温度、圧力による物性値変化を考慮すると共に、コ
ンピュータ等の情報処理装置を用いて高精度でかつ容易
に算出することを可能にした超臨界流体におけるピスト
ン効果による熱移動の算出方法を提供することにある。In view of the above-mentioned conventional problems, an object of the present invention is to consider changes in physical properties due to temperature and pressure, and to calculate with high accuracy and easily using an information processing device such as a computer. Another object of the present invention is to provide a method for calculating heat transfer due to the piston effect in a supercritical fluid described above.

【０００７】[0007]

【課題を解決するための手段】本発明は、上記の課題を
解決するために、次のような手段を採用した。The present invention employs the following means in order to solve the above problems.

【０００８】第１の手段は、超臨界流体におけるピスト
ン効果による熱移動の算出方法において、算出に必要な
解析領域、解析時間、初期条件、および境界条件の各情
報を入力する第１のステップと、初期温度Ｔ^０および初
期圧力ｐ^０に基づいて物性値を取得する第２のステップ
と、前記第２のステップまたは後述する第６のステップ
で取得された物性値を用いて、∂Ｔ／∂ｔ＝∇（λ
（Ｔ、ｐ）・∇Ｔ）／ρ（Ｔ、ｐ）・ｃ_ｐ（Ｔ、ｐ）式
に基づいて熱伝導による温度Ｔを算出する第３のステッ
プと、前記第２のステップまたは後述する第６のステッ
プで取得された物性値を用いて、∂ｐ／∂ｔ＝∫_Ｖρ
（Ｔ、ｐ）・α_ｐ（Ｔ、ｐ）・［∂Ｔ／∂ｔ］・ｄＶ／
∫_Ｖρ（Ｔ、ｐ）・χ_Ｔ（Ｔ、ｐ）・ｄＶ式に基づいて
圧力ｐを算出する第４のステップと、前記第２のステッ
プまたは後述する第６のステップで取得された物性値を
用いて、Ｑ＝［１−ｃ_Ｖ（Ｔ、ｐ）／ｃ_ｐ（Ｔ、ｐ）］
・［∂Ｔ／∂ｐ］_ρ（Ｔ、ｐ）・∂ｐ／∂ｔ式に基づい
て発熱量Ｑを算出する第５のステップと、前記第３のス
テップで算出された温度Ｔおよび前記第４のステップで
算出された圧力ｐに基づいて物性値を取得する第６のス
テップと、第５のステップで算出された発熱量Ｑおよび
第６のステップで取得された物性値を用いて、∂Ｔ／∂
ｔ＝∇（λ（Ｔ、ｐ）・∇Ｔ）／ρ（Ｔ、ｐ）・ｃ
_ｐ（Ｔ、ｐ）＋［１−ｃ _Ｖ（Ｔ、ｐ）／ｃ_ｐ（Ｔ、
ｐ）］・［∂Ｔ／∂ｐ］_ρ（Ｔ、ｐ）・∂ｐ／∂ｔ式に
基づいてピストン効果を考慮した温度Ｔを算出する第７
のステップとからなることを特徴とする。ただし、上記
の各物性値は、熱伝導率λ（Ｔ、ｐ）、定圧比熱ｃ
_ｐ（Ｔ、ｐ）、密度ρ（Ｔ、ｐ）、定積比熱ｃ_Ｖ（Ｔ、
ｐ）、体膨張率α_ｐ（Ｔ、ｐ）、等温圧縮率χ_Ｔ（Ｔ、
ｐ）、［∂Ｔ／∂ｐ］_ρ（Ｔ、ｐ）。The first means is a fix in a supercritical fluid.
In the method of calculating heat transfer due to
Analysis area, analysis time, initial condition, and boundary condition information
First step of inputting information and initial temperature T⁰And first
Period pressure p⁰Second step of obtaining physical property values based on
And the second step or the sixth step described later
∂T / ∂t = ∇ (λ
(T, p) ・ ∇T) / ρ (T, p) ・ c_p(T, p) formula
The third step for calculating the temperature T due to heat conduction based on
And a second step or a sixth step described later.
∂p / ∂t = ∫ using the physical property values acquired in_Vρ
(T, p) / α_p(T, p) ・ [∂T / ∂t] ・ dV /
∫_Vρ (T, p) · χ_TBased on (T, p) · dV formula
The fourth step of calculating the pressure p and the second step
Or the physical property values obtained in the sixth step described below
Using Q = [1-c_V(T, p) / c_p(T, p)]
・ [∂T / ∂p]_ρBased on (T, p) ・ ∂p / ∂t formula
The fifth step of calculating the heat generation amount Q by
In the temperature T calculated in step and the fourth step
A sixth step for obtaining a physical property value based on the calculated pressure p
Step and the calorific value Q calculated in the fifth step and
∂T / ∂ using the physical property values obtained in the sixth step
t = ∇ (λ (T, p) · ∇T) / ρ (T, p) · c
_p(T, p) + [1-c _V(T, p) / c_p(T,
p)] ・ [∂T / ∂p]_ρ(T, p) · ∂p / ∂t formula
7th which calculates the temperature T which considered the piston effect based on
And the steps of. However, above
The physical properties of are the thermal conductivity λ (T, p) and the constant pressure specific heat c
_p(T, p), density ρ (T, p), constant volume specific heat c_V(T,
p), body expansion rate α_p(T, p), isothermal compressibility χ_T(T,
p), [∂T / ∂p]_ρ(T, p).

【０００９】第２の手段は、第１の手段において、前記
第３のステップおよび前記第７のステップで算出された
温度Ｔについて、経時的に算出された最新の温度Ｔ
^ｎ＋１と経時的に前記最新の温度Ｔ^ｎ＋１の直前に算出
された温度Ｔ^ｎとに基づいて両者間の収束を判定し、収
束していないと判定されたときは前記第４のステップか
らの処理を繰り返させ、収束していると判定されたとき
は後述する第９のステップの処理を行わせる第８のステ
ップと、各ステップの算出結果を記憶ないし表示する第
９のステップとからなることを特徴とする。A second means is the latest temperature T calculated with time with respect to the temperature T calculated in the first step and the third step and the seventh step.
^{n + 1} and based on the time to the latest temperature T ^{n + 1} of the temperature T ⁿ calculated immediately before to determine the convergence between the two, when it is determined not converged processing from the fourth step Is repeated, and when it is determined that the convergence has occurred, an eighth step for performing the processing of the ninth step described later and a ninth step for storing or displaying the calculation result of each step. Characterize.

【００１０】第３の手段は、第２の手段において、前記
第９のステップの処理後、所定のタイムステップ数を処
理したかを判定し、処理していないと判定されたときは
前記第３のステップからの処理を行わせる第１０のステ
ップとからなることを特徴とする。The third means determines, in the second means, whether or not a predetermined number of time steps have been processed after the processing of the ninth step, and if it is determined that the predetermined number of time steps have not been processed, the third means. And a tenth step for performing the processing from step.

【００１１】第４の手段は、算出に必要な解析領域、解
析時間、初期条件、および境界条件の各情報を入力する
第１のステップと、初期温度Ｔ^０および初期圧力ｐ^０に
基づいて物性値を取得する第２のステップと、前記第２
のステップで取得された物性値を用いて、∂Ｔ／∂ｔ＝
∇（λ（Ｔ、ｐ）・∇Ｔ）／ρ（Ｔ、ｐ）・ｃ_ｐ（Ｔ、
ｐ）式に基づいて熱伝導による温度Ｔを算出する第３の
ステップと、前記第２のステップまたは後述する第６の
ステップで取得された物性値を用いて、∂ｐ／∂ｔ＝∫
_Ｖρ（Ｔ、ｐ）・α_ｐ（Ｔ、ｐ）・［∂Ｔ／∂ｔ］・ｄ
Ｖ／∫_Ｖρ（Ｔ、ｐ）・χ_Ｔ（Ｔ、ｐ）・ｄＶ式に基づ
いて圧力ｐを算出する第４のステップと、前記第２のス
テップまたは後述する第６のステップで取得された物性
値を用いて、Ｑ＝［１−ｃ_Ｖ（Ｔ、ｐ）／ｃ_ｐ（Ｔ、
ｐ）］・［∂Ｔ／∂ｐ］_ρ（Ｔ、ｐ）・∂ｐ／∂ｔ式に
基づいて発熱量Ｑを算出する第５のステップと、前記第
３のステップで算出された温度Ｔおよび前記第４のステ
ップで算出された圧力ｐに基づいて物性値を取得する第
６のステップと、第５のステップで算出された発熱量Ｑ
および第６のステップで取得された物性値を用いて、∂
Ｔ／∂ｔ＝∇（λ（Ｔ、ｐ）・∇Ｔ）／ρ（Ｔ、ｐ）・
ｃ_ｐ（Ｔ、ｐ）＋［１−ｃ_Ｖ（Ｔ、ｐ）／ｃ_ｐ（Ｔ、
ｐ）］・［∂Ｔ／∂ｐ］_ρ（Ｔ、ｐ）・∂ｐ／∂ｔ式に
基づいてピストン効果を考慮した温度Ｔを算出する第７
のステップと、前記第７のステップで算出された温度Ｔ
について、経時的に算出された最新の温度Ｔ^ｎ＋１と経
時的に前記最新の温度Ｔ^ｎ＋１の直前に算出された温度
Ｔ^ｎとに基づいて両者間の収束を判定し、収束していな
いと判定されたときは前記第４のステップからの処理を
繰り返させ、収束していると判定されたときは後述する
第９のステップの処理を行わせる第８のステップと、各
ステップの算出結果を記憶ないし表示する第９のステッ
プと、前記第９のステップの処理後、所定のタイムステ
ップ数を処理したかを判定し、処理していないと判定さ
れたときは前記第３のステップからの処理を行わせる第
１０のステップと、からなることを特徴とする超臨界流
体におけるピストン効果による熱移動を算出するプログ
ラム。ただし、上記の各物性値は、熱伝導率λ（Ｔ、
ｐ）、定圧比熱ｃ_ｐ（Ｔ、ｐ）、密度ρ（Ｔ、ｐ）、定
積比熱ｃ_Ｖ（Ｔ、ｐ）、体膨張率α_ｐ（Ｔ、ｐ）、等温
圧縮率χ_Ｔ（Ｔ、ｐ）、［∂Ｔ／∂ｐ］_ρ（Ｔ、ｐ）。A fourth means is to input the information of the analysis region, the analysis time, the initial condition and the boundary condition necessary for the calculation, and the physical properties based on the initial temperature T ⁰ and the initial pressure p ^0. A second step of obtaining a value, and the second step
∂T / ∂t = using the physical property values obtained in step
∇ (λ (T, p) · ∇T) / ρ (T, p) · c p (T,
p) ∂p / ∂t = ∫ using the third step of calculating the temperature T due to heat conduction based on the equation and the physical property values obtained in the second step or the sixth step described later.
_V ρ (T, p) ・ α _p (T, p) ・ [∂T / ∂t] ・ d
V / ∫ _V ρ (T, p) · χ _T (T, p) · dV obtained in the fourth step of calculating the pressure p based on the equation and the second step or the sixth step described later. by using the physical property _{values, Q = [1-c V} (T, p) / c p (T,
p)] · [∂T / ∂p] _ρ (T, p) · ∂p / ∂t based on the fifth step of calculating the calorific value Q, and the temperature T calculated in the third step. And a sixth step of obtaining a physical property value based on the pressure p calculated in the fourth step, and a calorific value Q calculated in the fifth step.
And using the physical property values obtained in the sixth step, ∂
T / ∂t = ∇ (λ (T, p) ・ ∇T) / ρ (T, p) ・
_{c p (T, p) +} [1-c V (T, p) / c p (T,
p)] · [∂T / ∂p] _ρ (T, p) · ∂p / ∂t based on the formula, the temperature T considering the piston effect is calculated.
And the temperature T calculated in the seventh step
For determines convergence between them on the basis of the temperature T ⁿ calculated immediately before the time the calculated latest temperature T ^{n + 1} and over time the latest temperature T ^{n + 1,} and does not converge determination When it is determined that the process from the fourth step is repeated, and when it is determined that the process is converged, an eighth step for performing the process of a ninth step described later and a calculation result of each step are stored. Or the ninth step of displaying and after the processing of the ninth step, it is determined whether or not a predetermined number of time steps have been processed, and when it is determined that the predetermined time step has not been processed, the processing from the third step is performed. A program for calculating heat transfer due to a piston effect in a supercritical fluid, characterized by comprising: a tenth step to be performed. However, the above-mentioned physical property values are the thermal conductivity λ (T,
p), pressure specific heat _c p (T, p), the density ρ (T, p), specific heat at constant volume _c V (T, p), volume expansion ratio α _p (T, _p), isothermal compressibility chi _T (T , P), [∂T / ∂p] _ρ (T, p).

【００１２】第５の手段は、コンピュータ読み取り可能
な記録媒体として、第４の手段に記載の超臨界流体にお
けるピストン効果による熱移動を算出するプログラムを
記録したことを特徴とする。A fifth means is characterized in that, as a computer-readable recording medium, a program for calculating the heat transfer due to the piston effect in the supercritical fluid described in the fourth means is recorded.

【００１３】[0013]

【発明の実施の形態】従来、Ｓｔｒａｕｂ等によって提
案され、超臨界流体におけるピストン効果による熱移動
現象を解析するための計算式としては下式のものが知ら
れている。DESCRIPTION OF THE PREFERRED EMBODIMENTS Conventionally, the following formula has been known as a calculation formula proposed by Straub et al. For analyzing a heat transfer phenomenon due to a piston effect in a supercritical fluid.

【００１４】 ∂Ｔ／∂ｔ＝∇（λ・∇Ｔ）／ρｃ_ｐ＋［１−ｃ_Ｖ／ｃ_ｐ］・ ［∂Ｔ／∂ｐ］_ρ・∂ｐ／∂ｔ・・・ （１）’ ∂ｐ／∂ｔ＝∫_Ｖρα_ｐ［∂Ｔ／∂ｔ］ｄＶ／∫_Ｖρχ_ＴｄＶ・・・ ・・（２）’ ここで、上記各式において、λは熱伝導率、ｃ_ｐは定圧
比熱、ρは密度、ｃ_Ｖは定積比熱（比熱比γ）、αは体
膨張率、ｃ_ｐは定圧比熱、χ_Ｔは等温圧縮率からなる物
性値である。[0014] ∂T / ∂t = ∇ (λ · ∇T) / ρc p + [1-c V / c p] · [∂T / ∂p] ρ · ∂p / ∂t ··· (1) '∂p / ∂t = ∫ _V ρα _p [∂T / ∂t] dV / ∫ _V ρχ _T dV ... (2)' where λ is thermal conductivity and c _p the specific heat at constant pressure, [rho is the density, c _V is specific heat at constant volume (specific heat ratio gamma), alpha is the body expansion, the c _p the specific heat at constant pressure, the chi _T is a physical property value consisting of isothermal compressibility.

【００１５】従来においても、これらの偏微分方程式を
適当な計算領域、境界条件、初期条件のもとに解くこと
により、超臨界流体におけるピストン効果による熱移動
現象を解明することが可能とされている。しかし、これ
らの偏微分方程式を解析的に解くのは一般的に困難なこ
とであり、また、本件発明者の知見によれば、超臨界流
体では温度、圧力の変化、特に、臨界点近傍では温度、
圧力のわずかな変化でも、これらの物性値に大きな変化
を与えることが分かり、各物性値が適正な値に設定され
なければ、実際の現象に即したピストン効果による熱移
動を算出することができない。Conventionally, it has been possible to solve the heat transfer phenomenon due to the piston effect in a supercritical fluid by solving these partial differential equations under appropriate calculation regions, boundary conditions and initial conditions. There is. However, it is generally difficult to solve these partial differential equations analytically, and according to the knowledge of the inventor of the present invention, in supercritical fluids, changes in temperature and pressure, especially in the vicinity of the critical point, temperature,
It has been found that even a slight change in pressure causes a large change in these physical property values, and unless each physical property value is set to an appropriate value, it is not possible to calculate the heat transfer due to the piston effect that matches the actual phenomenon. .

【００１６】そこで、本件発明においては、温度、圧力
の変化に伴って変化する物性値の変化を有効な手段によ
って逐次算出過程に反映させ、かつコンピュータ処理に
よって高精度でかつ容易に算出することを可能にした、
超臨界流体におけるピストン効果による熱移動の算出方
法を提案するものである。以下に、本発明の一実施形態
を図１ないし図６を用いて説明する。本発明では、上述
のごとく、実際の超臨界流体におけるピストン効果によ
る熱移動現象においては、各物性値が温度Ｔ、圧力ｐに
よって大きな影響を受けることから、各物性値を温度Ｔ
および圧力ｐの関数として表す。従って、各物性値であ
る、熱伝導率λはλ（Ｔ、ｐ）、定圧比熱ｃ_ｐはｃ
_ｐ（Ｔ、ｐ）、密度ρはρ（Ｔ、ｐ）、比熱比γはγ
（Ｔ、ｐ）（定積比熱ｃ_Ｖはｃ _Ｖ（Ｔ、ｐ））、体膨張
率αはα_ｐ（Ｔ、ｐ）、等温圧縮率χ_Ｔはχ_Ｔ（Ｔ、
ｐ）、そして、［∂Ｔ／∂ｐ］_ρは［∂Ｔ／∂ｐ］
_ρ（Ｔ、ｐ）として表わされ、上式（１）’式および
（２）’は下式のように表される。Therefore, in the present invention, temperature, pressure
The change of physical property value that changes with the change of
Is reflected in the sequential calculation process, and computer processing is performed.
Therefore, it is possible to calculate with high accuracy and easily.
Calculation method of heat transfer due to piston effect in supercritical fluid
It proposes a law. Below, one embodiment of the present invention
Will be described with reference to FIGS. In the present invention, the above
As a result, due to the piston effect in an actual supercritical fluid,
In the heat transfer phenomenon, each physical property value changes to temperature T and pressure p.
Therefore, the physical properties are greatly affected.
And as a function of pressure p. Therefore, in each physical property value
The thermal conductivity λ is λ (T, p) and the constant pressure specific heat c_pIs c
_p(T, p), density ρ is ρ (T, p), specific heat ratio γ is γ
(T, p) (constant volume specific heat c_VIs c _V(T, p)), body expansion
The rate α is α_p(T, p), isothermal compressibility χ_TIs χ_T(T,
p), and [∂T / ∂p]_ρIs [∂T / ∂p]
_ρIs expressed as (T, p) and is expressed by the above equation (1) ′ and
(2) 'is expressed by the following equation.

【００１７】 ∂Ｔ／∂ｔ＝∇（λ（Ｔ、ｐ）・∇Ｔ）／ρ（Ｔ、ｐ）・ｃ_ｐ（Ｔ、ｐ）＋ ［１−ｃ_Ｖ（Ｔ、ｐ）／ｃ_ｐ（Ｔ、ｐ）］・ ［∂Ｔ／∂ｐ］_ρ（Ｔ、ｐ）・∂ｐ／∂ｔ・・・（１） ∂ｐ／∂ｔ＝∫_Ｖρ（Ｔ、ｐ）・α_ｐ（Ｔ、ｐ）［∂Ｔ／∂ｔ］ｄＶ／ ∫_Ｖρ（Ｔ、ｐ）・χ_Ｔ（Ｔ、ｐ）・ｄＶ・・・（２）。[0017] ∂T / ∂t = ∇ (λ ( T, p) · ∇T) / ρ (T, p) · c p (T, p) + [1-c V (T, p) / c p (T, p)] ・ [∂T / ∂p] _ρ (T, p) ・ ∂p / ∂t ... (1) ∂p / ∂t = ∫ _V ρ (T, p) ・ α _p ( T, p) [∂T / ∂t ] dV / ∫ V ρ (T, p) · χ T (T, p) · dV ··· (2).

【００１８】以下に、上記の（１）式および（２）式に
基づいて、本件発明に係る超臨界流体におけるピストン
効果による熱移動の算出方法を図１を用いて説明する。A method for calculating heat transfer due to the piston effect in the supercritical fluid according to the present invention will be described below with reference to the above equations (1) and (2) with reference to FIG.

【００１９】図１は、本件発明に係る超臨界流体におけ
るピストン効果による熱移動の算出方法を示すフローチ
ャートである。FIG. 1 is a flow chart showing a method for calculating heat transfer due to the piston effect in a supercritical fluid according to the present invention.

【００２０】まず、ステップ１において、算出に必要な
解析領域、解析時間、初期条件、境界条件等の各種の情
報を設定し入力する。本実施形態では、超臨界流体とし
ては水を用いた。なお、一般に用いられる水は温度６４
７．０９６Ｋ、圧力２２．０６４ＭＰａで臨界点に達し
超臨界流体となる。算出には、図２に示すように、１次
元の物理モデルを採用した。First, in step 1, various kinds of information such as analysis area, analysis time, initial condition and boundary condition necessary for calculation are set and input. In this embodiment, water is used as the supercritical fluid. In addition, generally used water has a temperature of 64.
At 7.096 K and a pressure of 22.064 MPa, it reaches a critical point and becomes a supercritical fluid. For the calculation, as shown in FIG. 2, a one-dimensional physical model was adopted.

【００２１】計算領域は、同図に示すように、１次元で
長さは０．０１ｍであり、計算領域は１５３個のコント
ロールボリューム（なお、境界における厚さ０のコント
ロールボリュームも含む）で分割され、各コントロール
ボリュームの中心点に１５３点の格子点を設定した。な
お、本実施形態の物理モデルでは、一端付近の温度変化
が大きいので、この一端付近に多くの格子点を設けるた
めに、一端から長さ０．００１ｍの範囲に格子点の半数
を置いた。下式は、ｎ番目のコントロールボリュームの
厚さを表す式である。ｄ＝０．００１×（１．０−ｅ
^−α）ｅ^−ｎα／（１．０−ｅ^−７７α）ここで、αは
格子点をどのくらいの密度で配置するかの係数であり、
本実施形態では、α＝０．１とし、前記０．００１ｍの
位置から他端の０．０１ｍまではコントロールボリュー
ムを均等に配置した。As shown in the figure, the calculation area is one-dimensional and has a length of 0.01 m, and the calculation area is divided by 153 control volumes (including a control volume having a thickness of 0 at the boundary). Then, 153 grid points were set at the center points of the respective control volumes. In the physical model of the present embodiment, since the temperature change near one end is large, half of the grid points are placed within a range of 0.001 m from one end in order to provide many grid points near this one end. The following formula is a formula expressing the thickness of the n-th control volume. d = 0.001 * (1.0-e
^−α ) e ^−nα /(1.0−e ^−77α ) where α is a coefficient of how densely the lattice points are arranged,
In this embodiment, α = 0.1, and the control volume is evenly arranged from the position of 0.001 m to the other end of 0.01 m.

【００２２】問題としては、時間ｔ＝０において、前記
一端の温度が急に上昇し、定常状態において、領域全体
が前記一端の温度と同じになるまでのピストン効果によ
る熱移動を算出するものとし、境界条件としては、同図
に示すように、一端が温度一定（Ｔｈ）、他端が断熱
（∂Ｔ／∂ｘ＝０）という条件に設定し、初期条件は、
初期圧力が２２．１ＭＰａ、前記一端における初期温度
が６５１Ｋ、その他の領域は６５０Ｋとし、解析時間と
して、タイムステップは、０．１秒とし、タイムステッ
プ数は任意数設定した。The problem is that the heat transfer due to the piston effect is calculated until the temperature at the one end suddenly rises at time t = 0 and the entire region becomes the same as the temperature at the one end in the steady state. As the boundary condition, as shown in the figure, one end is set to a constant temperature (Th) and the other end is adiabatic (∂T / ∂x = 0), and the initial condition is
The initial pressure was 22.1 MPa, the initial temperature at the one end was 651 K, and the other region was 650 K. As the analysis time, the time step was set to 0.1 seconds, and the number of time steps was set to an arbitrary number.

【００２３】次に、ステップ２において、初期温度と初
期圧力から必要とされる物性値を取得する。ここで、初
期温度をＴ^０、初期圧力ｐ^０とすると、上記の各物性値
は、熱伝導率λ（Ｔ^０、ｐ^０）、定圧比熱ｃ_ｐ（Ｔ^０、
ｐ^０）、密度ρ（Ｔ^０、ｐ^０）、比熱比γ（Ｔ^０、
ｐ^０）（定積比熱ｃ_Ｖ（Ｔ^０、ｐ^０））、体膨張率α_ｐ
（Ｔ^０、ｐ^０）、等温圧縮率χ_Ｔ（Ｔ^０、ｐ^０）、［∂
Ｔ／∂ｐ］_ρ（Ｔ^０、ｐ^０）で表せる。なお、本実施形
態では、これらの物性値は、１９９７に作成されたＩＡ
ＰＷＳ実用国際状態式に基づいて作成された日本機械学
会蒸気表を参照して取得するように設定した。Next, in step 2, required physical property values are acquired from the initial temperature and the initial pressure. Here, the initial temperature ^{T 0,} when the initial pressure ^{p 0,} the physical property values described above, the thermal conductivity ^{λ (T 0, p 0)} , the constant-pressure specific heat _c p ^(T 0,
p ⁰ ), density ρ (T ⁰ , p ⁰ ), specific heat ratio γ (T ⁰ ,
p ⁰ ) (constant volume specific heat c _V (T ⁰ , p ⁰ )), body expansion coefficient α _p
(T ⁰ , p ⁰ ), isothermal compressibility χ _T (T ⁰ , p ⁰ ), [∂
T / ∂p] _ρ (T ⁰ , p ⁰ ). In addition, in the present embodiment, these physical property values are IA created in 1997.
It was set so as to be acquired by referring to the steam table of the Japan Society of Mechanical Engineers created based on the PWS practical international state equation.

【００２４】次に、ステップ３において、上記（１）式
において、前段の熱伝導項と後段の発熱項のうち、前段
の熱伝導項、∂Ｔ／∂ｔ＝∇（λ（Ｔ、ｐ）・∇Ｔ）／
ρ（Ｔ、ｐ）・ｃ_ｐ（Ｔ、ｐ）に基づいて、領域全体の
温度Ｔ^１を求める。ここで、物性値として、熱伝導率λ
（Ｔ^０、ｐ^０）、密度ρ（Ｔ^０、ｐ^０）、定圧比熱ｃ_ｐ
（Ｔ^０、ｐ^０）を用いると、上式は下式のように表され
る。 ∂Ｔ／∂ｔ＝∇（λ（Ｔ^０、ｐ^０）・∇Ｔ）／ρ（Ｔ^０、ｐ^０）・ ｃ_ｐ（Ｔ^０、ｐ^０）・・・（３−１） 温度Ｔ^１の算出に当たって、図３に示すように、コント
ロールボリュームの任意のコントロールボリュームにお
いて、任意の厚さ△ｘのコントロールボリュームおよび
その格子点Ｓ、およびそのコントロールボリュームに隣
接するコントロールボリュームの格子点Ｗ、Ｅ、格子点
Ｗと格子点Ｓ間および格子点Ｓと格子点Ｅ間の差分を
（δｘ）_ｗ、（δｘ）_ｅ、各コントロールボリューム間
の界面をｗおよびｅとし、格子点Ｓにおける過去（既
知）の温度、例えば、初期温度をＴ_Ｓ ^０とし、格子点
Ｗ、Ｓ、Ｅにおける温度をそれぞれ、Ｔ_Ｗ ^１、Ｔ_Ｓ ^１、
Ｔ_Ｅ ^１とすると、（３−１）式を離散化して表すと、下
式のように表される。 ａ_ＳＴ_Ｓ ^１＝ａ_ＥＴ_Ｅ ^１＋ａ_ＷＴ_Ｗ ^１＋ｂ・・・（４−１） ここで、ａ_Ｅ＝λ_ｅ（Ｔ^０、ｐ^０）／（δｘ）_ｅ、ａ_Ｗ
＝λ_ｗ（Ｔ^０、ｐ^０）／（δｘ）_ｗ、ａ_Ｓ ^０＝ρ
（Ｔ^０、ｐ^０）・ｃ_ｐ（Ｔ^０、ｐ^０）・△ｘ／△ｔ、ｂ
＝ａ_Ｓ・Ｔ_Ｓ ^０、ａ_Ｓ＝ａ_Ｅ＋ａ_Ｗ＋ａ_Ｓ ^０ 、△ｔは
タイムステップ間隔であり、密度ρ（Ｔ^０、ｐ^０）およ
び定圧比熱ｃ_ｐ（Ｔ^０、ｐ^０）はコントロールボリュー
ムＳの物性値を使用する。また熱伝導率λ_ｅ（Ｔ^０、ｐ
^０）、λ_ｗ（Ｔ ^０、ｐ^０）については、コントロールボ
リュームでの界面ｗ、ｅでの値が必要であり、実際に求
められているのは、コントロールボリュームＷ、Ｅの値
であるので、補間計数をｆ_ｅ≡（δｘ）_ｅ＋／（δｘ）
_ｅと定義すると、図４に示すように、界面ｅにおける熱
伝導率λ_ｅ（Ｔ^０、ｐ^０）は、λ_ｅ（Ｔ^０、ｐ^０）＝
（（１−ｆ_ｅ）／λ_Ｓ（Ｔ^０、ｐ^０）＋ｆ_ｅ／λ_Ｅ（Ｔ
^０、ｐ^０））で表される。従って、（４−１）式を利用
して、詳細は省略するが、領域全体の各格子点の温度を
未知数として連立方程式（マトリックス）を作成し、さ
らに、これに境界条件を入れて連立方程式を解くと、初
期温度Ｔ^０から熱伝導のみを考慮した、最初のタイムス
テップにおける領域全体の温度Ｔ^１が求められる。Next, in step 3, the above equation (1) is used.
Of the heat conduction term of the first stage and the heat generation term of the second stage,
Heat conduction term of ∂T / ∂t = ∇ (λ (T, p) ・ ∇T) /
ρ (T, p) · c_pBased on (T, p)
Temperature T¹Ask for. Here, as the physical property value, the thermal conductivity λ
(T⁰, P⁰), Density ρ (T⁰, P⁰), Constant pressure specific heat c_p
(T⁰, P⁰), The above equation is expressed as
It   ∂T / ∂t = ∇ (λ (T⁰, P⁰) ・ ∇T) / ρ (T⁰, P⁰) ・               c_p(T⁰, P⁰) ・ ・ ・ (3-1) Temperature T¹As shown in Fig. 3, when calculating
For any control volume of the roll volume
And a control volume of arbitrary thickness Δx and
Next to the grid point S and its control volume
Grid points W, E, grid points of the control volume that touches
Differences between W and grid point S and between grid point S and grid point E
(Δx)_w, (Δx)_eBetween each control volume
Let w and e be the interfaces of the
Knowledge) temperature, for example, the initial temperature is T_S ⁰And grid points
The temperatures at W, S, and E are T_W ¹, T_S ¹,
T_E ¹Then, when the equation (3-1) is discretized and expressed,
It is expressed as an expression.   a_ST_S ¹= A_ET_E ¹+ A_WT_W ¹+ B ... (4-1) Where a_E= Λ_e(T⁰, P⁰) / (Δx)_e, A_W
= Λ_w(T⁰, P⁰) / (Δx)_w, A_S ⁰= Ρ
(T⁰, P⁰) ・ C_p(T⁰, P⁰) ・ Δx / Δt, b
= A_S・ T_S ⁰, A_S= A_E+ A_W+ A_S ⁰ , Δt
Is the time step interval, and the density ρ (T⁰, P⁰) And
And constant pressure specific heat c_p(T⁰, P⁰) Is the control volume
Use the physical property value of S. And the thermal conductivity λ_e(T⁰, P
⁰), Λ_w(T ⁰, P⁰) For control
The values at the interfaces w and e in the ryme are necessary, and are actually calculated.
What is being set is the value of the control volume W, E
Therefore, the interpolation count is f_e≡ (δx)_{e +}/ (Δx)
_e, The heat at the interface e is defined as shown in FIG.
Conductivity λ_e(T⁰, P⁰) Is λ_e(T⁰, P⁰) =
((1-f_e) / Λ_S(T⁰, P⁰) + F_e/ Λ_E(T
⁰, P⁰)). Therefore, use equation (4-1)
Then, although details are omitted, the temperature of each grid point in the entire area is
Create simultaneous equations (matrix) as unknowns,
In addition, if the simultaneous equations are solved by adding boundary conditions to this,
Period temperature T⁰From the first time
Temperature T of the whole area at the step¹Is required.

【００２５】次に、ステップ４にて、上記（２）式、∂
ｐ／∂ｔ＝∫_Ｖρ（Ｔ、ｐ）・α_ｐ（Ｔ、ｐ）・［∂Ｔ
／∂ｔ］・ｄＶ／∫_Ｖρ（Ｔ、ｐ）・χ_Ｔ（Ｔ、ｐ）・
ｄＶを用いて圧力ｐを求める。（２）式は、物性値とし
て、各コントロールボリュームおける密度ρ（Ｔ^０、ｐ
^０）、体膨張率α_ｐ（Ｔ^０、ｐ^０）、等温圧縮率χ
_Ｔ（Ｔ^０、ｐ^０）、および各格子点における熱伝導のみ
を考慮した温度Ｔ^１に基づく［∂Ｔ^１／∂ｔ］を用いる
と、下式のように表される。 ∂ｐ／∂ｔ＝∫_Ｖρ（Ｔ^０、ｐ^０）α_ｐ（Ｔ^０、ｐ^０）・［∂Ｔ^１／∂ｔ］ ・ｄＶ／∫_Ｖρ（Ｔ^０、ｐ^０）・χ_Ｔ（Ｔ^０、ｐ^０）・ｄＶ・・ ・・（５−１） この（５−１）式を数値積分して圧力ｐ^１を求める。Next, in step 4, the above equation (2), ∂
_{p / ∂t = ∫ V ρ (} T, p) · α p (T, p) · [∂T
_{/ ∂t] · dV / ∫ V} ρ (T, p) · χ T (T, p) ·
The pressure p is determined using dV. Equation (2) is the physical property value of the density ρ (T ⁰ , p in each control volume).
⁰ ), body expansion rate α _p (T ⁰ , p ⁰ ), isothermal compression rate χ
_{Using T} (T ⁰ , p ⁰ ) and [∂T ¹ / ∂t] based on the temperature T ¹ in which only the heat conduction at each lattice point is taken into consideration, the following expression is obtained. ^{_{∂p / ∂t = ∫ V ρ (}} T 0, p 0) α p (T 0, p 0) · [∂T 1 / ∂t] · dV / ∫ V ρ (T 0, p 0) · χ T (T ⁰ , p ⁰ ) · dV ···· (5-1) This equation (5-1) is numerically integrated to obtain the pressure p ¹ .

【００２６】この式で特に注目すべき点は、圧力の時間
変化であり、領域内で温度が上がったことにより、圧力
がどれだけ上昇したかを把握するものである。ここで数
値積分は、ピストン効果を精度良く反映させるため、精
度の良いものでなければならない。また、温度が求めら
れている格子点の間隔が等間隔ではないことから、不等
間隔を扱えなければならない。従って、この計算におい
て不等間隔を扱うことのできるシンプソン（ｓｉｍｐｓ
ｏｎ）の公式を用いて、領域全体を積分して圧力を求め
る。The point of particular interest in this equation is the change in pressure over time, and it is to grasp how much the pressure has risen due to the rise in temperature within the region. Here, the numerical integration must be accurate because it accurately reflects the piston effect. In addition, since the lattice points for which the temperature is required are not evenly spaced, it is necessary to handle unequal spacing. Therefore, in this calculation, simpson (simps) that can handle unequal intervals
on) to calculate the pressure by integrating over the entire area.

【００２７】次に、ステップ５にて、上記（１）式の後
段の発熱項、Ｑ＝［１−ｃ_Ｖ／ｃ_ｐ］・［∂Ｔ／∂ｐ］
_ρ・∂ｐ／∂ｔから発熱量Ｑを求める。ここで、物性値
として、定積比熱ｃ_Ｖ（Ｔ^０、ｐ^０）、定圧比熱ｃ
_ｐ（Ｔ^０、ｐ^０）、密度ρ（Ｔ^０、ｐ^０）、比熱比γ
（Ｔ^０、ｐ^０）（定積比熱ｃ_Ｖ（Ｔ^０、ｐ^０））、体膨
張率α_ｐ（Ｔ^０、ｐ^０）、等温圧縮率χ_Ｔ（Ｔ^０、
ｐ^０）、［∂Ｔ／∂ｐ］_ρ（Ｔ^０、ｐ^０）を用いると、
発熱量Ｑ^１は、下式で表される。 Ｑ^１＝［１−ｃ_Ｖ（Ｔ^０、ｐ^０）／ｃ_ｐ（Ｔ^０、ｐ^０）］・ ［∂Ｔ／∂ｐ］_ρ（Ｔ^０、ｐ^０）・∂ｐ／∂ｔ・・・（６−１） ここで、圧力ｐは重力による影響の部分を除き、領域全
体で一定とし、本実施形態でも、圧力ｐは領域全体で一
定としている。Next, in step 5, above (1) in the subsequent stage of the heating _{section, Q = [1-c V} / c p] · [∂T / ∂p]
_The calorific value Q is calculated from _ρ · ∂p / ∂t. Here, as the physical property values, constant volume specific heat c _V (T ⁰ , p ⁰ ), constant pressure specific heat c
_p (T ⁰ , p ⁰ ), density ρ (T ⁰ , p ⁰ ), specific heat ratio γ
(T ⁰ , p ⁰ ) (constant volume specific heat c _V (T ⁰ , p ⁰ )), body expansion rate α _p (T ⁰ , p ⁰ ), isothermal compression rate χ _T (T ⁰ ,
p ⁰ ), [∂T / ∂p] _ρ (T ⁰ , p ⁰ )
The calorific value Q ¹ is represented by the following formula. ^{_{Q 1 = [1-c V}} (T 0, p 0) / c p (T 0, p 0)] · [∂T / ∂p] ρ (T 0, p 0) · ∂p / ∂t ·· (6-1) Here, the pressure p is constant over the entire region except for the portion affected by gravity, and in this embodiment also, the pressure p is constant over the entire region.

【００２８】次に、ステップ６では、ステップ３および
ステップ４にてそれぞれ求められた温度Ｔ^１および圧力
ｐ^１を用いて、各物性値を求める。これらの物性値はス
テップ２において求めたと同様に、ＩＡＰＷＳ実用国際
状態式に基づいて作成された日本機械学会蒸気表から、
熱伝導率λ（Ｔ^１、ｐ^１）、定圧比熱ｃ_ｐ（Ｔ^１、
ｐ^１）、密度ρ（Ｔ^１、ｐ^１）、比熱比γ（Ｔ^１、
ｐ^１）（定積比熱ｃ_Ｖ（Ｔ^１、ｐ^１）、体膨張率α
_ｐ（Ｔ^１、ｐ^１）、等温圧縮率χ _Ｔ（Ｔ^１、ｐ^１）、
［∂Ｔ／∂ｐ］_ρ（Ｔ^１、ｐ^１）を取得する。Next, in step 6, step 3 and
Temperature T obtained in step 4¹And pressure
p¹Each physical property value is calculated using. These physical properties are
Similar to the one obtained in Step 2, IAPWS Practical International
From the Japan Society of Mechanical Engineers steam table created based on the state equation,
Thermal conductivity λ (T¹, P¹), Constant pressure specific heat c_p(T¹,
p¹), Density ρ (T¹, P¹), Specific heat ratio γ (T¹,
p¹) (Constant volume specific heat c_V(T¹, P¹), Body expansion rate α
_p(T¹, P¹), Isothermal compressibility χ _T(T¹, P¹),
[∂T / ∂p]_ρ(T¹, P¹) To get.

【００２９】次に、ステップ７においてピストン効果を
考慮した温度Ｔを求める。熱伝導のみを考慮して求めら
れた温度Ｔ^１に対して、ピストン効果をも考慮した温度
をＴ^２とし、（１）式から温度Ｔ^２を求める。算出法
は、ステップ３で求めたと同様の手法によるが、ステッ
プ３で求めた方法と比べて、発熱量Ｑ^１が加算される点
で相違する。（１）式を離散化した式は、下式のように
表される。 ａ_ＳＴ_Ｓ ^２＝ａ_ＥＴ_Ｅ ^２＋ａ_ＷＴ_Ｗ ^２＋ｂ・・・（７−１） ここで、ａ_Ｅ＝λ_ｅ（Ｔ^１、ｐ^１）／（δｘ）_ｅ、ａ_Ｗ
＝λ_ｗ（Ｔ^１、ｐ^１）／（δｘ）_ｗ、ａ_Ｓ ^１＝ρ
（Ｔ^１、ｐ^１）・ｃ_ｐ（Ｔ^１、ｐ^１）・△ｘ／△ｔ、ｂ
＝ａ_Ｓ ^１・Ｔ_Ｓ ^１＋Ｑ^１Δｘ、ａ_Ｓ＝ａ_Ｅ＋ａ_Ｗ＋ａ_Ｓ
^１ 、△ｔはタイムステップ間隔である。ここでも詳細
は省略するが、（７−１）式を利用して、領域全体につ
いて連立方程式（マトリックス）を作成し、さらに、こ
れに境界条件を入れて連立方程式を解くと、温度Ｔ^１か
らピストン効果をも考慮した温度Ｔ^２が求められる。Next, in step 7, the temperature T considering the piston effect is obtained. With respect to the temperature T ¹ obtained by considering only the heat conduction, the temperature T ² which also takes the piston effect into consideration is set as T ^2, and the temperature T ² is obtained from the equation (1). The calculation method is the same as that obtained in step 3, but differs from the method obtained in step 3 in that the calorific value Q ¹ is added. The equation obtained by discretizing the equation (1) is expressed as the following equation. a _S T _S ² = a _E T _E ² + a _W T _W ² + b (7-1) Here, a _E = λ _e (T ¹ , p ¹ ) / (δx) _e , a _W
= Λ _w (T ¹ , p ¹ ) / (δx) _w , a _S ¹ = ρ
_{^{(T 1, p 1) ·}} c p (T 1, p 1) · △ x / △ t, b
= A _S ¹ · T _S ¹ + Q ¹ Δx, a _S = a _E + a _W + a _S
¹ and Δt are time step intervals. Although details are omitted here, if equations (7-1) are used to create simultaneous equations (matrix) for the entire region, and boundary equations are added to these to solve the simultaneous equations, the temperature T ¹ is calculated. The temperature T ² that takes into consideration the piston effect is also required.

【００３０】ピストン効果は、領域内のどこかでわずか
でも温度変化（上昇、下降）のあることが必要である
が、本発明の算出方法では熱伝導とピストン効果とによ
る熱移動現象を想定しており、この熱伝導による温度変
化がピストン効果を引き起こし、発熱項Ｑが領域全体に
影響を及ぼすものであるからピストン効果を表している
といえる。The piston effect requires a slight temperature change (rise, fall) anywhere in the region, but the calculation method of the present invention assumes a heat transfer phenomenon due to heat conduction and the piston effect. Therefore, it can be said that the temperature change due to the heat conduction causes the piston effect, and the heat generation term Q affects the entire region, so that it represents the piston effect.

【００３１】次に、ステップ８にて、ｎ＋１回目、ｎ回
目の計算でそれぞれ求められた温度Ｔ^ｎ＋１、Ｔ^ｎとす
るとき、所定値εに対して全ての格子点で｜（Ｔ^ｎ＋１
−Ｔ^ｎ）／Ｔ^ｎ｜＜εを判定する。温度Ｔ^２と温度Ｔ^１
間では収束することはあり得ないので、温度Ｔ^２と温度
Ｔ^１しか求められていない段階では、収束判定を行うこ
となく、ステップ４からの処理を行うようにしてもよ
い。その場合は、ステップ８での収束判定は｜（Ｔ^３−
Ｔ^２）／Ｔ^２｜＜εからとなる。なお、収束していると
判定されたときはステップ９の処理を行い、収束してい
ると判定されないときはステップ４からの処理を繰り返
す。なお、εとしては０．０１％（１０^−５）を用い
た。Next, in step 8, assuming that the temperatures T ^{n + 1} and T ⁿ are respectively obtained by the n + 1-th and n-th calculations, | (T ^{n + 1)} at all lattice points with respect to a predetermined value ε.
-T ⁿ ) / T ⁿ | <ε is determined. Temperature T ² and temperature T ¹
Since there is no possibility of convergence during the period, the processing from step 4 may be performed without performing the convergence determination at the stage where only the temperatures T ² and T ¹ are obtained. In that case, the convergence determination in step 8 is | (T ³ −
T ² ) / T ² | <ε. If it is determined that the convergence has occurred, the process of step 9 is performed, and if it is not determined that the convergence has occurred, the processes from step 4 are repeated. As ε, 0.01% (10 ⁻⁵ ) was used.

【００３２】ここで、収束していないと判定されたとき
のステップ４での処理は、先のステップ４での処理と同
様に、各コントロールボリュームおける密度ρ（Ｔ^１、
ｐ^１）、体膨張率α_ｐ（Ｔ^１、ｐ^１）、等温圧縮率χ_Ｔ
（Ｔ^１、ｐ^１）、および各格子点における温度Ｔ^２に基
づく［∂Ｔ^２／∂ｔ］を用い、下式を数値積分して、圧
力ｐ^２を求める。 ∂ｐ／∂ｔ＝∫_Ｖρ（Ｔ^１、ｐ^１）α_ｐ（Ｔ^１、ｐ^１）・［∂Ｔ^２／∂ｔ］・ ｄＶ／∫_Ｖρ（Ｔ^１、ｐ^１）・χ_Ｔ（Ｔ^１、ｐ^１）・ｄＶ・・ ・・（５−２）。Here, the processing in step 4 when it is determined that they have not converged is similar to the processing in step 4 above, the density ρ (T ¹ ,
p ¹ ), body expansion rate α _p (T ¹ , p ¹ ), isothermal compression rate χ _T
Using (T ¹ , p ¹ ) and [∂T ² / ∂t] based on the temperature T ² at each grid point, the following formula is numerically integrated to obtain the pressure p ² . ∂p / ∂t = ∫ _V ρ (T ¹ , p ¹ ) α _p (T ¹ , p ¹ ) ・ [∂T ² / ∂t] ・ dV / ∫ _V ρ (T ¹ , p ¹ ) ・ χ _T (T ¹ , p ¹ ) · dV ···· (5-2).

【００３３】次いで、ステップ５にて、先と同様にし
て、下式から発熱量Ｑ^２を求める。Ｑ^２＝［１−ｃ
_Ｖ（Ｔ^１、ｐ^１）／ｃ_ｐ（Ｔ^１、ｐ^１）］・ ［∂Ｔ／∂ｐ］_ρ（Ｔ^１、ｐ^１）・∂ｐ／∂ｔ・・・（６−２）。Then, in step 5, the calorific value Q ² is calculated from the following equation in the same manner as above. Q ² = [1-c
_{^{V (T 1, p 1)}} / c p (T 1, p 1)] · [∂T / ∂p] ρ (T 1, p 1) · ∂p / ∂t ··· (6-2).

【００３４】次いで、ステップ６にて、先と同様にし
て、ステップ３およびステップ４にてそれぞれ求められ
た温度Ｔ^２および圧力ｐ^２を用いて、各物性値を求め
る。これらの物性値はステップ２において求めたと同様
の手法により、熱伝導率λ（Ｔ^２、ｐ^２）、定圧比熱ｃ
_ｐ（Ｔ^２、ｐ^２）、密度ρ（Ｔ^２、ｐ^２）、比熱比γ
（Ｔ^２、ｐ^２）（定積比熱ｃ_Ｖ（Ｔ^２、ｐ^２）、体膨張
率α_ｐ（Ｔ^２、ｐ^２）、等温圧縮率χ_Ｔ（Ｔ^２、
ｐ^２）、［∂Ｔ／∂ｐ］_ρ（Ｔ^２、ｐ^２）を求める。Then, in step 6, each physical property value is obtained using the temperature T ² and the pressure p ² obtained in step 3 and step 4, respectively, in the same manner as above. These physical property values are calculated by the same method as that obtained in step 2, the thermal conductivity λ (T ² , p ² ), the constant pressure specific heat c
_p (T ² , p ² ), density ρ (T ² , p ² ), specific heat ratio γ
(T ² , p ² ) (constant volume specific heat c _V (T ² , p ² ), body expansion rate α _p (T ² , p ² ), isothermal compression rate χ _T (T ² ,
p ² ) and [∂T / ∂p] _ρ (T ² , p ² ) are obtained.

【００３５】さらに、ステップ７にて、先と同様に、下
式からピストン効果を考慮した温度Ｔ^３を求める。 ａ_ＳＴ_Ｓ ^３＝ａ_ＥＴ_Ｅ ^３＋ａ_ＷＴ_Ｗ ^３＋ｂ・・・（７−２） ここで、ａ_Ｅ＝λ_ｅ（Ｔ^２、ｐ^２）／（δｘ）_ｅ、ａ_Ｗ
＝λ_ｗ（Ｔ^２、ｐ^２）／（δｘ）_ｗ、ａ_Ｓ ^２＝ρ
（Ｔ^２、ｐ^２）・ｃ_ｐ（Ｔ^２、ｐ^２）・△ｘ／△ｔ、ｂ
＝ａ_Ｓ ^２Ｔ_Ｓ ^２＋Ｑ^２Δｘ、ａ_Ｓ＝ａ_Ｅ＋ａ_Ｗ＋ａ_Ｓ ^２
となり、詳細は省略するが、（７−２）式を利用し
て、領域全体について連立方程式（マトリックス）を作
成し、さらに、これに境界条件を入れて連立方程式を解
くと、温度Ｔ^２からピストン効果をも考慮した温度Ｔ^３
が求められる。Further, in step 7, similarly to the above, the temperature T ³ considering the piston effect is obtained from the following equation. a _S T _S ³ = a _E T _E ³ + a _W T _W ³ + b (7-2) where a _E = λ _e (T ² , p ² ) / (δx) _e , a _W
= Λ _w (T ² , p ² ) / (δx) _w , a _S ² = ρ
_{^{(T 2, p 2) ·}} c p (T 2, p 2) · △ x / △ t, b
= A _S ² T _S ² + Q ² Δx, a _S = a _E + a _W + a _S ²
Although the details are omitted, if equations (7-2) are used to create simultaneous equations (matrix) for the entire region, and then the simultaneous equations are solved by including boundary conditions in these, the temperature T ² is calculated. Temperature T ³ that also considers the piston effect
Is required.

【００３６】次いで、ステップ８で再び、所定値εに対
して、全ての格子点で｜（Ｔ^３−Ｔ^２）／Ｔ^２｜＜εを
判定することにより、温度Ｔ^３と温度Ｔ^２間の収束を判
定する。収束していないと判定されたときは、上記に述
べたように、再度ステップ４からの処理を繰り返すが、
このような処理を繰り返して、ステップ８での収束判定
において、全ての格子点で｜（Ｔ^ｎ＋１−Ｔ^ｎ）／Ｔ^ｎ
｜＜εを満足していると判定された時は、ステップ９に
移行する。Then, in step 8, again, for a predetermined value ε, by determining | (T ³ −T ² ) / T ² | <ε at all the lattice points, the temperature between the temperatures T ³ and T ² is determined. Determine the convergence of. If it is determined that they have not converged, as described above, the process from step 4 is repeated,
By repeating such processing, in the convergence determination in step 8, | (T ^{n + 1} −T ⁿ ) / T ⁿ at all grid points
If it is determined that | <ε is satisfied, the process proceeds to step 9.

【００３７】ステップ９では、各ステップ、例えば、ス
テップ８での処理の結果算出された各コントロールボリ
ュームにおける全ての格子点の温度Ｔ^ｎ＋１を記憶し、
また必要に応じて、任意のタイムステップ毎にグラフ表
示する。In step 9, the temperatures T ^{n + 1} of all the lattice points in each control volume calculated as a result of the processing in each step, for example, step 8 are stored,
If necessary, a graph is displayed for each arbitrary time step.

【００３８】次に、ステップ１０では、次のタイムステ
ップの処理を行うか否かを判定し、次のタイムステップ
での処理を行う場合は、それまで算出された算出結果を
初期値として、ステップ３において、再び熱伝導のみの
新たな温度Ｔ^１を算出するための処理を開始する。設定
された所定のタイムステップ数の処理を終えたとき全て
の処理を終了する。Next, at step 10, it is determined whether or not the processing at the next time step is to be performed. When the processing at the next time step is to be performed, the calculation result calculated up to that point is used as an initial value to perform the step. In 3, the process for calculating the new temperature T ¹ only for heat conduction is started again. When the processing of the set predetermined number of time steps is completed, all the processing is ended.

【００３９】図５は、本発明の超臨界流体のピストン効
果による熱移動の算出方法によって算出された結果を示
す図である。同図において、横軸は領域（０．０００〜
０．０１０ｍ）の位置ｘ（ｍ）、縦軸は温度Ｔ（Ｋ）、
ｔ（秒）＝０．１，１．０，１０，１００はそれぞれ算
出処理の経過時間を表す。FIG. 5 is a diagram showing the results calculated by the method for calculating the heat transfer by the piston effect of the supercritical fluid of the present invention. In the figure, the horizontal axis represents the region (0.000 to 0).
0.010 m) position x (m), the vertical axis represents temperature T (K),
t (seconds) = 0.1, 1.0, 10, 100 respectively represent the elapsed time of the calculation process.

【００４０】これは、超臨界水（臨界点温度Ｔｃ＝６４
７．０９６Ｋ、圧力ｐｃ＝２２．０６４ＭＰａ）では、
初期温度Ｔ^０＝６５０Ｋ、初期圧力ｐ^０＝２２．１ＭＰ
ａ、ｔ＝０において一端（ｘ＝０．０００）の温度が急
に６５１Ｋになった場合の温度分布の時間変化を示すも
のであり、ピストン効果は熱の伝わりが速く、領域全体
が同じ温度になるという特性があるが、この結果は、そ
の状況を良く表している。特に、ｔ＝０．１秒後でも、
ｔ＝１．０秒後でも他端（ｘ＝０．０１０）で温度上昇
が見られ、これは熱伝導だけでは考えられない現象であ
り、また、境界から徐々に熱が伝わるのではなくて、領
域全体に温度が伝わるという傾向が見られ、これは微小
重力実験で得られた結果とよく一致し、本発明の超臨界
流体の熱移動の算出方法が、ピストン効果による熱移動
の極めて有効な算出方法であることが解る。This is due to supercritical water (critical point temperature Tc = 64
At 7.096 K and pressure pc = 22.064 MPa),
Initial temperature T ⁰ = 650K, initial pressure p ⁰ = 22.1MP
It shows the change over time in the temperature distribution when the temperature at one end (x = 0.000) suddenly reaches 651K at a and t = 0. The piston effect shows that heat transfer is fast and the entire region has the same temperature. However, this result shows the situation well. In particular, even after t = 0.1 seconds,
Even after t = 1.0 seconds, a temperature rise was observed at the other end (x = 0.010), which is a phenomenon that cannot be considered only by heat conduction, and heat is not gradually transferred from the boundary. , There is a tendency that the temperature is transmitted to the entire region, which is in good agreement with the result obtained in the microgravity experiment, and the calculation method of the heat transfer of the supercritical fluid of the present invention is extremely effective for the heat transfer by the piston effect. It turns out that this is a simple calculation method.

【００４１】なお、図６は本発明の超臨界流体のピスト
ン効果による熱移動の算出方法を実現するためのハード
ウェア構成の一例を示す図である。同図において、１は
コンピュータプログラムに従って各種の演算を実行する
中央演算処理装置、２は演算処理に必要なデータ、例え
ば、熱移動の算出に必要な解析領域、解析時間、初期条
件、境界条件等の各種の情報を入力する入力装置、３
は、超臨界流体のピストン効果による熱移動の算出を実
行するプログラム等を記憶するメインメモリ、４は超臨
界流体のピストン効果による熱移動の算出において必要
なデータ、例えば、算出された温度Ｔや圧力ｐに基づい
て物性値を取得するために参照される日本機械学会蒸気
表等のデータや、超臨界流体の熱移動算出の結果が記憶
されるファイル装置、５は超臨界流体のピストン効果に
よる熱移動の算出結果等を画像表示したりプリントアウ
トする出力装置である。FIG. 6 is a diagram showing an example of a hardware configuration for realizing the method for calculating the heat transfer by the piston effect of the supercritical fluid according to the present invention. In the figure, 1 is a central processing unit for executing various calculations according to a computer program, 2 is data required for the calculation processing, for example, an analysis region required for calculation of heat transfer, analysis time, initial conditions, boundary conditions, etc. Input device for inputting various information of 3
Is a main memory for storing a program or the like for calculating the heat transfer due to the piston effect of the supercritical fluid, and 4 is data necessary for calculating the heat transfer due to the piston effect of the supercritical fluid, such as the calculated temperature T or File device 5 in which data such as the steam table of the Japan Society of Mechanical Engineers referred to for obtaining the physical property value based on the pressure p and the result of heat transfer calculation of the supercritical fluid are stored This is an output device that displays an image or prints out the calculation result of heat transfer and the like.

【００４２】[0042]

【発明の効果】請求項１ないし請求項３に記載の発明に
よれば、はじめに熱伝導を算出し、それに基づいてピス
トン効果を算出するようにしているので、ピストン効果
が、わずかでも温度上昇（ないし温度下降）しなければ
発生しない性質を見極めた算出方法が得られる。また、
温度、圧力の変化に伴って変化する物性値の変化を逐次
算出過程に反映させ、かつコンピュータ処理によって高
精度でかつ容易に算出するようにしたので、工業化レベ
ルでの超臨界流体におけるピストン効果による熱移動の
算出が容易となり、各種の工業分野での応用が可能とな
る。According to the invention described in claims 1 to 3, the heat transfer is first calculated and the piston effect is calculated based on the calculated heat transfer. It is possible to obtain a calculation method that identifies the properties that do not occur unless the temperature decreases). Also,
Changes in physical properties that change with changes in temperature and pressure are sequentially reflected in the calculation process, and are calculated with high precision and ease by computer processing. Calculation of heat transfer becomes easy, and application in various industrial fields becomes possible.

【００４３】請求項２に記載の発明によれば、収束判定
を行いながら、ピストン効果による熱移動を算出してい
るので、正確な算出結果が得られる。According to the second aspect of the present invention, since the heat transfer due to the piston effect is calculated while performing the convergence determination, an accurate calculation result can be obtained.

【００４４】請求項３に記載の発明によれば、タイムス
テップ数を適正な値に設定することにより、各時間毎の
詳細な超臨界流体におけるピストン効果による熱移動に
よる詳細な温度分布を得ることができる。According to the third aspect of the present invention, by setting the number of time steps to an appropriate value, it is possible to obtain a detailed temperature distribution by heat transfer due to the piston effect in the detailed supercritical fluid at each time. You can

【００４５】請求項４記載の発明によれば、はじめに熱
伝導を算出し、それに基づいてピストン効果を算出する
ようにしているので、ピストン効果が、わずかでも温度
上昇（ないし温度下降）しなければ発生しない性質を見
極めた算出プログラムが得られる。また、温度、圧力の
変化に伴って変化する物性値の変化を逐次算出過程に反
映させ、また収束判定を行いながら、ピストン効果によ
る熱移動を算出しているので、ピストン効果による熱移
動の正確な算出結果が得られる。その結果、工業化レベ
ルでの超臨界流体におけるピストン効果による熱移動の
算出が容易となり、各種の工業分野での応用が可能とな
る。According to the fourth aspect of the present invention, the heat conduction is first calculated and the piston effect is calculated based on the heat conduction. Therefore, if the piston effect does not increase in temperature at all (or decrease in temperature). It is possible to obtain a calculation program that identifies properties that do not occur. In addition, changes in physical properties that change with changes in temperature and pressure are sequentially reflected in the calculation process, and heat transfer due to the piston effect is calculated while performing convergence determination, so the heat transfer due to the piston effect is accurate. It is possible to obtain various calculation results. As a result, it becomes easy to calculate the heat transfer due to the piston effect in the supercritical fluid at the industrial level, and the application in various industrial fields becomes possible.

【００４６】請求項５記載の発明によれば、記録媒体に
本発明に係るプログラムを記録して頒布することによ
り、広く本発明の利用が可能となる。According to the invention of claim 5, the program according to the present invention is recorded in a recording medium and distributed, whereby the present invention can be widely used.

【図面の簡単な説明】[Brief description of drawings]

【図１】本件発明に係る超臨界流体におけるピストン効
果による熱移動の算出方法を示すフローチャートであ
る。FIG. 1 is a flowchart showing a method for calculating heat transfer due to a piston effect in a supercritical fluid according to the present invention.

【図２】格子点、コントロールボリューム、タイムステ
ップを示した超臨界流体のピストン効果による熱移動を
算出するための物理モデルを示す図である。FIG. 2 is a diagram showing a physical model for calculating heat transfer due to a piston effect of a supercritical fluid showing lattice points, control volumes, and time steps.

【図３】格子点Ｗ、Ｓ、Ｅとコントロールボリュームと
の関係を示す図である。FIG. 3 is a diagram showing a relationship between grid points W, S, E and a control volume.

【図４】格子点Ｓ、Ｅとコントロールボリュームの界面
ｅとの関係を示す図である。FIG. 4 is a diagram showing a relationship between lattice points S and E and an interface e of a control volume.

【図５】本発明に係る超臨界流体のピストン効果による
熱移動の算出方法によって算出された温度分布の時間変
化を示す図である。FIG. 5 is a diagram showing a time change of a temperature distribution calculated by a method of calculating heat transfer due to a piston effect of a supercritical fluid according to the present invention.

【図６】本発明の超臨界流体におけるピストン効果によ
る熱移動の算出方法を実現するためのハードウェア構成
の一例を示す図である。FIG. 6 is a diagram showing an example of a hardware configuration for realizing a method of calculating heat transfer due to a piston effect in a supercritical fluid of the present invention.

【符号の説明】[Explanation of symbols]

１ 中央演算処理装置 ２ 入力装置 ３ メインメモリ ４ ファイル装置 ５ 出力装置 1 Central processing unit 2 input devices 3 main memory 4 file device 5 Output device

───────────────────────────────────────────────────── フロントページの続き (72)発明者 相澤 崇史 宮城県仙台市宮城野区苦竹４−２−１ 独 立行政法人産業技術総合研究所 東北セン ター内 (72)発明者 金久保 光央 宮城県仙台市宮城野区苦竹４−２−１ 独 立行政法人産業技術総合研究所 東北セン ター内 (72)発明者 斎藤 功夫 宮城県仙台市宮城野区苦竹４−２−１ 独 立行政法人産業技術総合研究所 東北セン ター内 Ｆターム(参考） 2G040 AA01 AB08 AB12 BA22 HA15 HA16 5B056 BB00 HH00    ─────────────────────────────────────────────────── ─── Continued front page    (72) Inventor Takashi Aizawa             4-2-1 Maitake, Miyagino-ku, Sendai City, Miyagi Prefecture, Germany             National Institute of Advanced Industrial Science and Technology Tohoku Sen             Inside (72) Inventor Mitsuo Kanakubo             4-2-1 Maitake, Miyagino-ku, Sendai City, Miyagi Prefecture, Germany             National Institute of Advanced Industrial Science and Technology Tohoku Sen             Inside (72) Inventor Isao Saito             4-2-1 Maitake, Miyagino-ku, Sendai City, Miyagi Prefecture, Germany             National Institute of Advanced Industrial Science and Technology Tohoku Sen             Inside F-term (reference) 2G040 AA01 AB08 AB12 BA22 HA15                       HA16                 5B056 BB00 HH00

Claims (5)

【特許請求の範囲】[Claims] 【請求項１】 算出に必要な解析領域、解析時間、初期
条件、および境界条件の各情報を入力する第１のステッ
プと、 初期温度Ｔ^０および初期圧力ｐ^０に基づいて物性値を取
得する第２のステップと、 前記第２のステップで取得された物性値を用いて、 ∂Ｔ／∂ｔ＝∇（λ（Ｔ、ｐ）・∇Ｔ）／ρ（Ｔ、ｐ）
・ｃ_ｐ（Ｔ、ｐ） に基づいて熱伝導による温度Ｔを算出する第３のステッ
プと、 前記第２のステップまたは後述する第６のステップで取
得された物性値を用いて、 ∂ｐ／∂ｔ＝∫_Ｖρ（Ｔ、ｐ）・α_ｐ（Ｔ、ｐ）・［∂
Ｔ／∂ｔ］・ｄＶ／∫_Ｖρ（Ｔ、ｐ）・χ_Ｔ（Ｔ、ｐ）
・ｄＶ に基づいて圧力ｐを算出する第４のステップと、 前記第２のステップまたは後述する第６のステップで取
得された物性値を用いて、 Ｑ＝［１−ｃ_Ｖ（Ｔ、ｐ）／ｃ_ｐ（Ｔ、ｐ）］・［∂Ｔ
／∂ｐ］_ρ（Ｔ、ｐ）・∂ｐ／∂ｔ に基づいて発熱量Ｑを算出する第５のステップと、 前記第３のステップで算出された温度Ｔおよび前記第４
のステップで算出された圧力ｐに基づいて物性値を取得
する第６のステップと、 第５のステップで算出された発熱量Ｑおよび第６のステ
ップで取得された物性値を用いて、 ∂Ｔ／∂ｔ＝∇（λ（Ｔ、ｐ）・∇Ｔ）／ρ（Ｔ、ｐ）
・ｃ_ｐ（Ｔ、ｐ）＋［１−ｃ_Ｖ（Ｔ、ｐ）／ｃ_ｐ（Ｔ、
ｐ）］・［∂Ｔ／∂ｐ］_ρ（Ｔ、ｐ）・∂ｐ／∂ｔ に基づいてピストン効果を考慮した温度Ｔを算出する第
７のステップと、 からなることを特徴とする超臨界流体におけるピストン
効果による熱移動の算出方法。 ただし、上記の各物性値は、熱伝導率λ（Ｔ、ｐ）、定
圧比熱ｃ_ｐ（Ｔ、ｐ）、密度ρ（Ｔ、ｐ）、定積比熱ｃ
_Ｖ（Ｔ、ｐ）、体膨張率α_ｐ（Ｔ、ｐ）、等温圧縮率χ
_Ｔ（Ｔ、ｐ）、［∂Ｔ／∂ｐ］_ρ（Ｔ、ｐ）。1. A physical property value is obtained based on a first step of inputting each information of an analysis region, an analysis time, an initial condition, and a boundary condition necessary for calculation, and an initial temperature T ⁰ and an initial pressure p ^0. Using the second step and the physical property values obtained in the second step, ∂T / ∂t = ∇ (λ (T, p) · ∇T) / ρ (T, p)
∂p / using the third step of calculating the temperature T due to heat conduction based on c _p (T, p) and the physical property values obtained in the second step or the sixth step described later. _{∂t = ∫ V ρ (T,} p) · α p (T, p) · [∂
T / ∂t] ・ dV / ∫ _V ρ (T, p) ・ χ _T (T, p)
- a fourth step of calculating the pressure p on the basis of the dV, using the sixth to physical property values obtained in the step of the second step or later, Q = [1-c V (T, p) _{/ c p (T, p)} ] · [∂T
/ ∂p] _ρ (T, p) · ∂p / ∂t based on the fifth step, the temperature T calculated in the third step and the fourth step
∂T using the sixth step of acquiring the physical property value based on the pressure p calculated in the step, the calorific value Q calculated in the fifth step, and the physical property value acquired in the sixth step. / ∂t = ∇ (λ (T, p) ・ ∇T) / ρ (T, p)
_{· C p (T, p)} + [1-c V (T, p) / c p (T,
p)] · [∂T / ∂p] _ρ (T, p) · ∂p / ∂t based on the seventh step of calculating the temperature T considering the piston effect, and Calculation method of heat transfer due to piston effect in critical fluids. However, the physical properties of the above, thermal conductivity lambda (T, p), pressure specific heat _c p (T, p), the density ρ (T, p), specific heat at constant volume c
_V (T, p), body expansion rate α _p (T, p), isothermal compression rate χ
_T (T, p), [∂T / ∂p] _ρ (T, p). 【請求項２】 前記第７のステップで算出された温度Ｔ
について、経時的に算出された最新の温度Ｔ^ｎ＋１と経
時的に前記最新の温度Ｔ^ｎ＋１の直前に算出された温度
Ｔ^ｎとに基づいて両者間の収束を判定し、収束していな
いと判定されたときは前記第４のステップからの処理を
繰り返させ、収束していると判定されたときは後述する
第９のステップの処理を行わせる第８のステップと、 各ステップの算出結果を記憶ないし表示する第９のステ
ップと、 からなることを特徴とする請求項１に記載の超臨界流体
の熱移動算出方法。2. The temperature T calculated in the seventh step
For determines convergence between them on the basis of the temperature T ⁿ calculated immediately before the time the calculated latest temperature T ^{n + 1} and over time the latest temperature T ^{n + 1,} and does not converge determination When it is determined that the process from the fourth step is repeated, and when it is determined that the process is converged, an eighth step for performing the process of a ninth step described later and a calculation result of each step are stored. Or a ninth step of displaying, and the method for calculating heat transfer of a supercritical fluid according to claim 1, wherein: 【請求項３】 前記第９のステップの処理後、所定のタ
イムステップ数を処理したかを判定し、処理していない
と判定されたときは前記第３のステップからの処理を行
わせる第１０のステップと、 からなることを特徴とする請求項２に記載の超臨界流体
におけるピストン効果による熱移動の算出方法。3. After the processing of the ninth step, it is determined whether or not a predetermined number of time steps have been processed, and when it is determined that the predetermined number of time steps have not been processed, the processing from the third step is performed. The method for calculating heat transfer due to the piston effect in a supercritical fluid according to claim 2, characterized in that: 【請求項４】 算出に必要な解析領域、解析時間、初期
条件、および境界条件の各情報を入力する第１のステッ
プと、 初期温度Ｔ^０および初期圧力ｐ^０に基づいて物性値を取
得する第２のステップと、 前記第２のステップで取得された物性値を用いて、 ∂Ｔ／∂ｔ＝∇（λ（Ｔ、ｐ）・∇Ｔ）／ρ（Ｔ、ｐ）
・ｃ_ｐ（Ｔ、ｐ） に基づいて熱伝導による温度Ｔを算出する第３のステッ
プと、 前記第２のステップまたは後述する第６のステップで取
得された物性値を用いて、 ∂ｐ／∂ｔ＝∫_Ｖρ（Ｔ、ｐ）・α_ｐ（Ｔ、ｐ）・［∂
Ｔ／∂ｔ］・ｄＶ／∫_Ｖρ（Ｔ、ｐ）・χ_Ｔ（Ｔ、ｐ）
・ｄＶ に基づいて圧力ｐを算出する第４のステップと、 前記第２のステップまたは後述する第６のステップで取
得された物性値を用いて、 Ｑ＝［１−ｃ_Ｖ（Ｔ、ｐ）／ｃ_ｐ（Ｔ、ｐ）］・［∂Ｔ
／∂ｐ］_ρ（Ｔ、ｐ）・∂ｐ／∂ｔ に基づいて発熱量Ｑを算出する第５のステップと、 前記第３のステップで算出された温度Ｔおよび前記第４
のステップで算出された圧力ｐに基づいて物性値を取得
する第６のステップと、 第５のステップで算出された発熱量Ｑおよび第６のステ
ップで取得された物性値を用いて、 ∂Ｔ／∂ｔ＝∇（λ（Ｔ、ｐ）・∇Ｔ）／ρ（Ｔ、ｐ）
・ｃ_ｐ（Ｔ、ｐ）＋［１−ｃ_Ｖ（Ｔ、ｐ）／ｃ_ｐ（Ｔ、
ｐ）］・［∂Ｔ／∂ｐ］_ρ（Ｔ、ｐ）・∂ｐ／∂ｔ に基づいてピストン効果を考慮した温度Ｔを算出する第
７のステップと、 前記第７のステップで算出された温度Ｔについて、経時
的に算出された最新の温度Ｔ^ｎ＋１と経時的に前記最新
の温度Ｔ^ｎ＋１の直前に算出された温度Ｔ^ｎとに基づい
て両者間の収束を判定し、収束していないと判定された
ときは前記第４のステップからの処理を繰り返させ、収
束していると判定されたときは後述する第９のステップ
の処理を行わせる第８のステップと、 各ステップの算出結果を記憶ないし表示する第９のステ
ップと、 前記第９のステップの処理後、所定のタイムステップ数
を処理したかを判定し、処理していないと判定されたと
きは前記第３のステップからの処理を行わせる第１０の
ステップと、 からなることを特徴とする超臨界流体におけるピストン
効果による熱移動を算出するプログラム。ただし、上記
の各物性値は、熱伝導率λ（Ｔ、ｐ）、定圧比熱ｃ
_ｐ（Ｔ、ｐ）、密度ρ（Ｔ、ｐ）、定積比熱ｃ_Ｖ（Ｔ、
ｐ）、体膨張率α_ｐ（Ｔ、ｐ）、等温圧縮率χ_Ｔ（Ｔ、
ｐ）、［∂Ｔ／∂ｐ］_ρ（Ｔ、ｐ）。4. A physical property value is acquired based on a first step of inputting each information of an analysis region, an analysis time, an initial condition, and a boundary condition necessary for calculation, and an initial temperature T ⁰ and an initial pressure p ^0. Using the second step and the physical property values obtained in the second step, ∂T / ∂t = ∇ (λ (T, p) · ∇T) / ρ (T, p)
∂p / using the third step of calculating the temperature T due to heat conduction based on c _p (T, p) and the physical property values obtained in the second step or the sixth step described later. _{∂t = ∫ V ρ (T,} p) · α p (T, p) · [∂
T / ∂t] ・ dV / ∫ _V ρ (T, p) ・ χ _T (T, p)
- a fourth step of calculating the pressure p on the basis of the dV, using the sixth to physical property values obtained in the step of the second step or later, Q = [1-c V (T, p) _{/ c p (T, p)} ] · [∂T
/ ∂p] _ρ (T, p) · ∂p / ∂t based on the fifth step, the temperature T calculated in the third step and the fourth step
∂T using the sixth step of acquiring the physical property value based on the pressure p calculated in the step, the calorific value Q calculated in the fifth step, and the physical property value acquired in the sixth step. / ∂t = ∇ (λ (T, p) ・ ∇T) / ρ (T, p)
_{· C p (T, p)} + [1-c V (T, p) / c p (T,
p)] · [∂T / ∂p] _ρ (T, p) · ∂p / ∂t based on the seventh step of calculating the temperature T considering the piston effect, and the seventh step. and the temperature T, to determine the convergence between the two on the basis of the temperature T ⁿ calculated immediately before the time the calculated latest temperature T ^{n + 1} and over time the latest temperature T ^{n + 1,} have converged If it is determined that there is not, the process from the fourth step is repeated, and if it is determined that the process has converged, an eighth step that causes the process of the ninth step described later to be performed, and calculation of each step. A ninth step of storing or displaying the result, and after the processing of the ninth step, it is determined whether or not a predetermined number of time steps have been processed, and when it is determined that the processing has not been performed, the third step is executed. The 10th step When a program for calculating the heat transfer by the piston effect in the supercritical fluid, characterized in that it consists of. However, the above-mentioned physical property values are thermal conductivity λ (T, p), constant pressure specific heat c
_p (T, _p), the density ρ (T, p), specific heat at constant volume _c V (T,
p), the body expansion rate α _p (T, p), the isothermal compression rate χ _T (T,
p), [∂T / ∂p] _ρ (T, p). 【請求項５】 請求項４に記載のプログラムを記録した
コンピュータ読み取り可能な記録媒体。5. A computer-readable recording medium in which the program according to claim 4 is recorded.