Information processing method, physical property information providing system, and information processing program

Abstract

Provided is an information processing method which is performed by a computer so as to identify a gas-liquid interface of a gas-liquid mixed phase flow flowing in a flow path, the information processing method comprising: acquiring first gas-liquid mixed phase flow information including the free energy of each of one or more substances that are included in the gas-liquid mixed phase flow; acquiring reaction information related to a reaction occurring between one or more substances that are included in the gas phase and one or more substances that are included in the liquid phase; identifying a gas-liquid interface of the gas-liquid mixed phase flow by the free energy of the gas phase and the free energy of the liquid phase; and outputting gas-liquid interface information indicating the identified gas-liquid interface.

Description

Information processing method, physical property information providing system, and information processing program

[0001] The present disclosure relates to a technique for simulating a gas-liquid multiphase flow.

[0002] Patent Document 1 discloses an analysis device that calculates the flow velocity, pressure, and enthalpy of a gas-liquid two-phase fluid for each of a plurality of compartments based on the phase field method, calculates the mass flow fraction in each of a plurality of compartments including the gas-liquid interface based on the pressure and enthalpy, calculates physical property values ​​such as the density, viscosity, thermal conductivity, and specific heat of the two phases based on the pressure, enthalpy, mass flow fraction, and mathematical formula information, and outputs the calculation results.

[0003] However, the technology of Patent Document 1 does not capture changes in component concentration due to dissolution of the gas phase into the liquid phase and chemical reactions, so further improvement is needed to accurately predict the gas-liquid interface in gas-liquid multiphase flows.

[0004] JP 2012-247094 A

[0005] The present disclosure has been made to solve such problems, and aims to provide a technology for accurately predicting the gas-liquid interface in a gas-liquid multiphase flow.

[0006] An information processing method in one aspect of the present disclosure is an information processing method performed by a computer to identify a gas-liquid interface of a gas-liquid multiphase flow flowing in a flow path, and includes: acquiring first gas-liquid multiphase flow information including the free energy of each of one or more substances included in the gas-liquid multiphase flow; acquiring reaction information related to reactions occurring between one or more substances included in the gas phase and one or more substances included in the liquid phase; identifying the gas-liquid interface of the gas-liquid multiphase flow using the free energy of the gas phase and the free energy of the liquid phase; and outputting gas-liquid interface information indicating the identified gas-liquid interface.

[0007] According to the present disclosure, the gas-liquid interface of a gas-liquid multiphase flow can be accurately predicted.

[0008] 10 is a block diagram showing an example of the overall configuration of an information processing system in an embodiment. FIG. 11 is a block diagram showing a detailed configuration of a calculation process unit. FIG. 12 is a flowchart showing an example of processing of the information processing system in an embodiment. FIG. 13 is a flowchart showing an example of initial condition acquisition processing. FIG. 14 is a flowchart showing an example of pre-process processing. FIG. 15 is a flowchart showing details of calculation process calculation processing. FIG. 16 is a flowchart showing details of calculation process calculation processing continued from FIG. 6. FIG. 17 is a flowchart showing details of calculation process calculation processing continued from FIG. 7. FIG. 18 is a diagram showing an example of a flow channel model. FIG. 19 is a diagram of the flow channel model in which mesh points are set. FIG. 10 is a diagram showing an example of an input screen. FIG. 11 is a diagram showing an example of an output screen. FIG. 12 is a diagram showing an example of an output screen. FIG. 13 is a diagram showing simulation results of a comparative example. FIG. 14 is a diagram showing an example of a physical property list. FIG. 15 is a diagram showing another example of a physical property list.

[0009] (Findings underlying the present disclosure) The phase field method, which has a thickness of the gas-liquid interface between the gas and liquid and smoothly changes physical properties, is known as one method for determining the gas-liquid interface of a gas-liquid multiphase flow. The phase field method is a method for determining phase change interfaces, such as solid-liquid interfaces and gas-liquid interfaces, based on the principle of free energy minimization. If the physical properties of each phase can be precisely described, it is possible in principle to precisely determine the gas-liquid interface. For example, as shown in Patent Document 1, the phase field method is used as a method for describing the gas-liquid interface in the boiling phenomenon of water.

[0010] However, the technology disclosed in Patent Document 1 does not capture the reaction property values ​​of chemical reactions, and therefore further improvement is required to accurately predict the gas-liquid interface of a gas-liquid multiphase flow.

[0011] Therefore, the present inventors have come up with the idea of ​​the present disclosure based on the finding that the gas-liquid interface can be calculated precisely by incorporating the parameters necessary for the chemical reaction in advance, calculating the reaction property values ​​of the gas-liquid multiphase flow based on the reaction information, and calculating the gas-liquid interface of the gas-liquid multiphase flow using the calculated reaction property values.

[0012] (1) An information processing method according to one aspect of the present disclosure is an information processing method performed by a computer to identify a gas-liquid interface of a gas-liquid multiphase flow flowing through a flow path, the information processing method including: acquiring first gas-liquid multiphase flow information including the free energy of each of one or more substances included in the gas-liquid multiphase flow; acquiring reaction information regarding reactions occurring between one or more substances included in the gas phase and one or more substances included in the liquid phase; identifying the gas-liquid interface of the gas-liquid multiphase flow using the free energy of the gas phase and the free energy of the liquid phase; and outputting gas-liquid interface information indicating the identified gas-liquid interface.

[0013] According to this configuration, first gas-liquid multiphase flow information including the free energy of each of one or more substances in the gas-liquid multiphase flow is acquired, and reaction information regarding a reaction occurring between one or more substances in the gas phase and one or more substances in the liquid phase is acquired. This allows parameters necessary for a chemical reaction to be acquired in advance. Then, the gas-liquid interface of the gas-liquid multiphase flow is identified based on the free energy of the gas phase and the free energy of the liquid phase. This configuration allows for accurate prediction of the gas-liquid interface.

[0014] (2) The information processing method described in (1) above may further include acquiring structural information including a structure including the flow path and boundary conditions for calculating the gas-liquid interface of the gas-liquid multiphase flow, and calculating one or more reaction property values ​​obtained by the reaction based on the first gas-liquid multiphase flow information, the structural information, and the reaction information.

[0015] According to this configuration, the gas-liquid interface of the gas-liquid multiphase flow is identified by further using structural information including the structure including the flow path and the boundary conditions for calculating the gas-liquid interface of the gas-liquid multiphase flow.

[0016] (3) In the information processing method described in (2) above, calculating the gas-liquid interface may include calculating the gas-liquid interface by executing a simulation using a phase field method with the structural information, the reaction information, and the one or more reaction property values ​​as input, and the simulation using the phase field method may include describing the free energy as a sum of many-body potentials based on a cluster expansion.

[0017] In this case, the gas-liquid interface can be predicted more precisely.

[0018] (4) In the information processing method described in any one of (1) to (3) above, the reaction may include dissolving the one or more substances contained in the gas phase into the one or more substances in the liquid phase, and identifying the gas-liquid interface may include calculating the gas-liquid interface so that the difference between the free energy of the gas phase and the free energy of the liquid phase is small.

[0019] In this case, the reaction includes dissolution, and the gas-liquid interface is calculated so that the difference between the free energy of the gas phase and the free energy of the liquid phase is small, so the gas-liquid interface can be predicted more accurately.

[0020] (5) In the information processing method described in any one of (2) to (4) above, identifying the gas-liquid interface may include calculating a surface tension based on the temperature or dissolution calculated based on the reaction information, and identifying the gas-liquid interface so that the surface tension matches the external pressure.

[0021] (6) In the information processing method described in any one of (2) to (5) above, calculating the one or more reaction property values ​​may further include calculating at least one of the amount of reaction component and the heat of reaction obtained by the reaction based on the first gas-liquid multiphase flow information, the structural information, and the reaction information, and calculating the temperature of the gas-liquid multiphase flow resulting from the reaction based on the first gas-liquid multiphase flow information, the structural information, and at least one of the amount of reaction component and the heat of reaction.

[0022] According to this configuration, at least one of the amount of reaction components and the heat of reaction, and the temperature of the gas-liquid multiphase flow are calculated as one or more reaction property values, so that the gas-liquid interface can be predicted more accurately by taking into account changes in components due to chemical reactions.

[0023] (7) In the information processing method described in any one of (2) to (6) above, the method may further include calculating second gas-liquid multiphase flow information including one or more second property values ​​of the gas-liquid multiphase flow based on the first gas-liquid multiphase flow information, the calculated one or more reactant property values, and the calculated gas-liquid interface, wherein calculating the one or more reactant property values ​​includes calculating the one or more reactant property values ​​obtained by the reaction based on the second gas-liquid multiphase flow information, the structural information, and the reaction information, and calculating the gas-liquid interface of the gas-liquid multiphase flow may include calculating the gas-liquid interface of the gas-liquid multiphase flow based on the second gas-liquid multiphase flow information, the structural information, the reaction information, and the one or more reactant property values.

[0024] According to this configuration, second gas-liquid multiphase flow information including one or more second property values ​​of the gas-liquid multiphase flow is calculated based on the first gas-liquid multiphase flow information, the calculated one or more reactant property values, and the calculated gas-liquid interface. Then, the calculated second gas-liquid multiphase flow information is used to further calculate one or more reactant property values. Furthermore, the calculated one or more reactant property values ​​are used to calculate the gas-liquid interface. This allows for more accurate prediction of the gas-liquid interface.

[0025] (8) In the information processing method described in any one of (1) to (7) above, the output gas-liquid interface information may indicate a gas-liquid interface calculated at a predetermined elapsed time.

[0026] According to this configuration, it is possible to output the gas-liquid interface information when a predetermined time has elapsed since the start of the simulation.

[0027] (9) In the information processing method described in any one of (2) to (8) above, acquiring the structural information may include acquiring a plurality of mesh points included in the structure, and may further include outputting second gas-liquid multiphase flow information including one or more second physical property values ​​of the gas-liquid multiphase flow calculated for each of the plurality of mesh points.

[0028] According to this configuration, the second gas-liquid multiphase flow information calculated for each of a plurality of mesh points can be output.

[0029] (10) In the information processing method described in any one of (1) to (9) above, the method may further include outputting second gas-liquid multiphase flow information including one or more second physical property values ​​of the gas-liquid multiphase flow calculated at a predetermined elapsed time in a simulation.

[0030] According to this configuration, it is possible to output the second gas-liquid multiphase flow information for each mesh point when a predetermined time has elapsed since the start of the simulation.

[0031] (11) In the information processing method described in any one of (2) to (10) above, identifying the gas-liquid interface may include calculating the first gas-liquid multiphase flow information, the structural information, the reaction information, and the one or more reaction property values ​​by substituting them into a Cahn-Hilliard or Allen-Cahn equation.

[0032] According to this configuration, the gas-liquid interface is calculated using the Cahn-Hilliard or Allen-Cahn equation, so that the gas-liquid interface can be predicted accurately.

[0033] (12) In the information processing method described in any one of (1) to (11) above, acquiring the structural information may include acquiring a plurality of mesh points included in the structure, and may further include determining whether or not a calculated value of the gas-liquid interface of the gas-liquid multiphase flow calculated for each of the plurality of mesh points is below a predetermined reference value, and changing a convergence parameter used in calculating the gas-liquid interface of the gas-liquid multiphase flow in accordance with the determination result, and identifying the gas-liquid interface may include identifying the gas-liquid interface of the gas-liquid multiphase flow based on the changed convergence parameter.

[0034] According to this configuration, the calculation of the gas-liquid interface is repeated while changing the convergence parameters used in the calculation of the gas-liquid interface until the calculated value of the gas-liquid interface falls below a predetermined reference value, thereby enabling the solution of the gas-liquid interface to converge quickly.

[0035] (13) In another aspect of the present disclosure, a physical property information providing system is a physical property information providing system for a gas-liquid multiphase flow flowing through a flow path, and includes: a display control unit that causes a display unit to display a first image including a flow path model representing the flow path and data representing mesh points present in the flow path model or time information indicating a predetermined elapsed time; and an image generation unit that, upon receiving a selection of the data representing the mesh points or the time information, generates a second image including data representing calculated values ​​of one or more physical properties of the gas-liquid multiphase flow.

[0036] According to this configuration, it is possible to display the calculated values ​​of one or more physical properties at any mesh point or any elapsed time.

[0037] (14) In the physical property information providing system described in (13) above, when the image generating unit receives a selection of data indicating the mesh point, it may generate the second image including data representing calculated values ​​of one or more physical properties of the gas-liquid multiphase flow at the selected mesh point.

[0038] According to this configuration, it is possible to display the calculated values ​​of one or more physical properties at any mesh point.

[0039] (15) In another aspect of the present disclosure, an information processing program is an information processing program for identifying a gas-liquid interface of a gas-liquid multiphase flow flowing in a flow path, and causes a computer to execute the following operations: acquire first gas-liquid multiphase flow information including the free energy of each of one or more substances included in the gas-liquid multiphase flow; acquire reaction information regarding reactions occurring between one or more substances included in the gas phase and one or more substances included in the liquid phase; identify the gas-liquid interface of the gas-liquid multiphase flow using the free energy of the gas phase and the free energy of the liquid phase; and output gas-liquid interface information indicating the identified gas-liquid interface.

[0040] According to this configuration, it is possible to provide an information processing program that can precisely calculate the gas-liquid interface.

[0041] (16) In another aspect of the present disclosure, an information processing system is an information processing system for calculating a gas-liquid interface of a gas-liquid multiphase flow flowing in a flow path, and includes: an acquisition unit that acquires first gas-liquid multiphase flow information including the free energy of each of one or more substances included in the gas-liquid multiphase flow; and reaction information regarding reactions occurring between one or more substances included in the gas phase and one or more substances included in the liquid phase; a calculation process unit that identifies the gas-liquid interface of the gas-liquid multiphase flow using the free energy of the gas phase and the free energy of the liquid phase; and an output unit that outputs gas-liquid interface information indicating the identified gas-liquid interface.

[0042] This configuration makes it possible to provide an information processing system that can precisely identify the gas-liquid interface.

[0043] The present disclosure can provide a physical property information providing program that causes a computer to function as a physical property information providing system, or a physical property information providing method having a characteristic configuration of a physical property information providing system. The present disclosure can also be realized as an information processing system operated by an information processing program. Furthermore, these physical property information providing programs and information processing programs can be distributed on computer-readable non-transitory recording media such as CD-ROMs or via communication networks such as the Internet.

[0044] Note that each of the embodiments described below represents a specific example of the present disclosure. The numerical values, shapes, components, steps, and step orders shown in the following embodiments are merely examples and are not intended to limit the present disclosure. Furthermore, among the components in the following embodiments, components that are not described in the independent claims that represent the highest concept are described as optional components. Furthermore, in all of the embodiments, the respective contents can be combined.

[0045] 1 is a block diagram showing an example of the overall configuration of an information processing system 1 according to an embodiment. The information processing system 1 is configured with a computer that calculates the gas-liquid interface of a gas-liquid multiphase flow flowing through a flow path. The information processing system 1 includes a processor 10, an operation unit 20, a display unit 30, and a memory 40.

[0046] The processor 10 is composed of a central processing unit (CPU). The processor 10 includes an initial condition acquisition unit 11, a calculation process unit 12, and an output unit 13. The output unit 13 includes a display control unit 14, an image generation unit 15, and a data transfer unit 16. The initial condition acquisition unit 11 to the output unit 13 are each realized by the CPU executing an information processing program stored in memory 40. However, this is just one example, and the initial condition acquisition unit 11 to the output unit 13 may be composed of dedicated hardware circuits. In the example of FIG. 1, the information processing system 1 is composed of a single computer, but it may also be composed of multiple computers. In this case, the initial condition acquisition unit 11 to the output unit 13 may be distributed across multiple computers. The information processing system 1 may also be composed of a cloud server. In this case, the cloud server may output processing results to a terminal (not shown).

[0047] The initial condition acquisition unit 11 acquires first gas-liquid multiphase flow information, structural information, and reaction information. The first gas-liquid multiphase flow information includes the free energy of each of one or more substances included in the gas-liquid multiphase flow. In detail, the first gas-liquid multiphase flow information includes one or more first physical property values ​​of the gas-liquid multiphase flow, the type of one or more substances included in the gas-liquid multiphase flow, and free energy information of each of the one or more substances.

[0048] The first physical property values ​​include, for example, an initial value of pressure, an initial value of flow velocity, an initial value of temperature, an initial value of concentration, an initial value of constant pressure heat capacity, an initial value of thermal conductivity, an initial value of viscosity coefficient, an initial value of density, an initial value of latent heat of vaporization, an initial value of gas-liquid interface, and an initial value of diffusion coefficient. In this embodiment, calculated values ​​are calculated at multiple mesh points set in the structure. Therefore, the first physical property values ​​include values ​​set at each mesh point.

[0049] The one or more types of substances are, for example, the names of substances that constitute the gas-liquid multiphase flow.

[0050] The free energy information includes an initial value of the free energy of each of one or more substances. Furthermore, if the reaction includes dissolution, the free energy information also includes an initial value of the solvation free energy. The initial values ​​of the free energy and the solvation free energy are set for each of the plurality of mesh points.

[0051] The structural information includes a flow path model (an example of a structure) including the flow path and boundary conditions for calculating the gas-liquid interface of the gas-liquid multiphase flow. The flow path model is a three-dimensional model of the flow path through which the gas-liquid multiphase flow flows.

[0052] FIG. 9 is a diagram showing an example of a flow channel model 1000. The flow channel model 1000 is composed of a tubular member having a ring-shaped cross section. The flow channel model 1000 includes a mechanism unit 1001 and a flow channel 1002. The mechanism unit 1001 is a substantial unit that defines the flow channel 1002. The flow channel 1002 is a space through which a gas-liquid multiphase flow flows. The flow channel model 1000 has an x-axis in the width direction, a y-axis in the height direction, and a z-axis in the longitudinal direction. The flow channel model 1000 includes an inlet 1005 and an outlet 1006. The gas-liquid multiphase flow flows through the flow channel 1002 from the inlet 1005 to the outlet 1006.

[0053] 10 is a view of the flow path model 1000 of FIG. 9 as viewed from the inlet 1005 side. The flow path model 1000 includes a heat source 1003. The heat source 1003 is disposed, for example, at a predetermined position on a central axis O1 in the longitudinal direction of the flow path 1002. The heat source 1003 heats the gas-liquid multiphase flow flowing through the flow path 1002.

[0054] The boundary conditions are determined by the inner wall of the mechanism part 1001, the inlet 1005, and the outlet 1006. Furthermore, the boundary conditions are also determined by the heat source 1003.

[0055] Fig. 11 is a diagram showing a flow channel model 1000 in which mesh points are set. The flow channel model 1000 is a three-dimensional model created in advance in a computer space. As shown in Fig. 11, a plurality of mesh points P are set in the flow channel model 1000. The mesh points P are set three-dimensionally in a flow channel 1002. The mesh points P are also set three-dimensionally inside a mechanical section 1001.

[0056] 15 is a cross-sectional view of the flow channel model 1000 as seen from the side. It can be seen that the heat source 1003 is installed on the central axis of the flow channel model 1000. By adding heat from the heat source 1003 to the gas-liquid multiphase flow, it is possible to simulate the behavior of the gas-liquid multiphase flow when it is heated, for example.

[0057] The reaction information is information about a reaction occurring between one or more substances contained in the gas phase and one or more substances contained in the liquid phase. The reaction information includes a chemical reaction formula, a reaction rate constant, molar enthalpy, an initial value of the latent heat of vaporization, and an initial value of the heat of reaction. The chemical reaction formula is a chemical reaction formula that indicates a chemical reaction occurring in a gas-liquid multiphase flow.

[0058] The calculation process unit 12 identifies the gas-liquid interface of the gas-liquid multiphase flow based on the free energy of the gas phase and the free energy of the liquid phase.

[0059] The calculation processing unit 12 calculates one or more reaction property values ​​obtained by the reaction based on the first gas-liquid multiphase flow information, the structure information, and the reaction information. The reaction property values ​​include the amounts of reaction components and the temperature. The calculation processing unit 12 identifies the gas-liquid interface of the gas-liquid multiphase flow based on the first gas-liquid multiphase flow information, the structure information, the reaction information, and the one or more reaction property values.

[0060] Specifically, the calculation unit 12 calculates the gas-liquid interface using a simulation based on the phase field method, which inputs structural information, reaction information, and one or more reaction property values. Here, the calculation unit 12 uses the phase field method, in which the gas-liquid interface is described by free energy described as the sum of many-body potentials based on cluster expansion. The process of calculating the gas-liquid interface is executed by the gas-liquid interface calculation unit 127, which will be described later.

[0061] The reaction involves dissolving one or more substances contained in the gas phase into one or more substances in the liquid phase.

[0062] The calculation process unit 12 specifies the gas-liquid interface so that the difference between the free energy of the gas phase and the free energy of the liquid phase becomes small.

[0063] The calculation process unit 12 calculates one or more reaction property values ​​by performing a process of calculating at least one of the amounts of reactants and the heat of reaction obtained by the reaction based on the first gas-liquid multiphase flow information, the structure information, and the reaction information, and a process of calculating the temperature of the gas-liquid multiphase flow resulting from the reaction based on the first gas-liquid multiphase flow information, the structure information, and at least one of the amounts of reactants and the heat of reaction. The process of calculating at least one of the amounts of reactants and the heat of reaction is performed by a chemical reaction calculation unit 123 described below. The process of calculating the temperature of the gas-liquid multiphase flow resulting from the reaction is performed by a heat transfer calculation unit 124 described below.

[0064] The calculation processing unit 12 calculates second gas-liquid multiphase flow information including one or more second property values ​​of the gas-liquid multiphase flow based on the first gas-liquid multiphase flow information, the calculated one or more reactant property values, and the identified gas-liquid interface. In this case, the calculation processing unit 12 calculates one or more reactant property values ​​obtained by the reaction based on the second gas-liquid multiphase flow information, the structure information, and the reaction information. Furthermore, the calculation processing unit 12 calculates the gas-liquid interface of the gas-liquid multiphase flow based on the second gas-liquid multiphase flow information, the structure information, the reaction information, and the one or more reactant property values. The process of calculating the second gas-liquid multiphase flow information is executed by the fluid calculation unit 122, the diffusion calculation unit 125, and the property calculation unit 126, which will be described later.

[0065] For example, the second physical property values ​​are the first physical property values ​​updated by the simulation, and thus, like the first physical property values, the second physical property values ​​include pressure, flow velocity, temperature, constant pressure heat capacity, thermal conductivity, latent heat of vaporization, diffusion coefficient, viscosity coefficient, and density.

[0066] The output unit 13 generates gas-liquid interface information indicating the gas-liquid interface calculated at a predetermined elapsed time, and outputs the generated gas-liquid interface information. For example, the output unit 13 may display an output screen 1600 ( FIG. 16 ) of the generated gas-liquid interface information on the display unit 30.

[0067] The initial condition acquisition unit 11 may acquire, as the structural information, a plurality of mesh points included in the structure. In this case, the output unit 13 may output second gas-liquid multiphase flow information including one or more second physical property values ​​of the gas-liquid multiphase flow calculated for each of the plurality of mesh points. For example, the output unit 13 may display an output screen 1300 ( FIG. 13 ) of the second gas-liquid multiphase flow information on the display unit 30. The plurality of mesh points are calculation points at which various values ​​of the simulation are calculated.

[0068] The output unit 13 may output second gas-liquid multiphase flow information including one or more second physical property values ​​of the gas-liquid multiphase flow calculated at a predetermined elapsed time by the calculation process unit 12. For example, the output unit 13 may display an output screen 1400 of the second gas-liquid multiphase flow information on the display unit 30.

[0069] The output unit 13 includes a display control unit 14 , an image generation unit 15 , and a data transfer unit 16 .

[0070] The display control unit 14 causes the display unit 30 to display a first image including the flow channel model 1000 included in the structure information and data representing mesh points present in the flow channel model 1000 or time information indicating a predetermined elapsed time. An example of the first image is an input screen 1200 ( FIG. 12 ).

[0071] FIG. 12 is a diagram showing an example of an input screen 1200. The input screen 1200 includes a flow path model 1000 and a time input field 1201. The input screen 1200 displays the flow path model 1000 as viewed from the side. The user can use the operation unit 20 to input an instruction to select an arbitrary position P1 in the flow path model 1000 displayed on the input screen 1200. Alternatively, the user can input an arbitrary elapsed time (time information) from the start to the end of the simulation in the time input field 1201 included in the input screen 1200. Note that although mesh points are not displayed in the flow path model 1000 displayed on the input screen 1200, mesh points are preset in the flow path model 1000 as shown in FIG. 11 .

[0072] When the image generation unit 15 receives a selection of data representing mesh points on the input screen 1200, it generates a second image including data representing calculated values ​​of one or more physical properties of the gas-liquid multiphase flow and displays the generated second image on the display unit 30. For example, when an instruction to select position P1 is input, the image generation unit 15 displays a second image including data representing calculated values ​​of one or more physical properties of the gas-liquid multiphase flow at position P1. An example of the second image displayed in this case is the output screen 1300 ( FIG. 13 ). In the example of the output screen 1300, the velocity of the gas-liquid multiphase flow is used as the calculated value of the physical property.

[0073] 13 is a diagram showing an example of an output screen 1300. In the output screen 1300, the vertical axis indicates the velocity corresponding to the arc length of a mesh point cross section including a selected mesh point, and the horizontal axis indicates the arc length in the selected mesh point cross section. The arc length refers to the length of an arc from a reference point on the circumference of the mesh point cross section. The mesh point cross section refers to a disk cut out from the cross section of the flow channel 1002 at the selected mesh point, and the disk has a radius equal to the distance between the mesh point and the central axis of the flow channel 1002 in the longitudinal direction. The reference point can be, for example, a vertex in the height direction in the mesh point cross section.

[0074] Five graphs corresponding to elapsed times of 0 s, 0.1 s, 0.2 s, 0.3 s, and 0.4 s are displayed on the output screen 1300. For example, when one focuses on the graph for 0 s, the distribution of speed at 0 s according to the arc length can be seen.

[0075] The output screen 1300 displays the velocity of the selected mesh point, but may also display the flow rate and concentration of the selected mesh point.

[0076] Note that the selection of position P1 on input screen 1200 is not limited to a mesh point. When a position between mesh points (sub-mesh point) is selected as position P1, data representing the calculated value at position P1 is calculated by interpolating data representing the calculated values ​​at the surrounding mesh points. Then, a second image including data representing the interpolated calculated value is displayed on display unit 30.

[0077] When the image generation unit 15 receives a selection of time information on the input screen 1200, it generates a second image including data representing calculated values ​​of one or more physical properties of the gas-liquid multiphase flow. For example, when a time is input in the time input field 1201, the image generation unit 15 generates a second image including data representing calculated values ​​at the input time and displays the generated second image on the display unit 30. An example of the second image displayed in this case is the output screen 1400 or the output screen 1600. Note that the time input in the time input field 1201 may be a numerical input. Alternatively, the time input field 1201 may be configured with a slide bar for selecting an arbitrary time step in the simulation period. In this case, the user selects an arbitrary time by sliding the slide bar.

[0078] If a sub-time step between two time steps is selected, the output unit 13 may interpolate the calculated value of the sub-time step using the data before and after the time step. Examples of the interpolation method that can be used include linear approximation, quadratic approximation, polynomial approximation, and asymptotic expansion.

[0079] FIG. 14 is a diagram showing an example of an output screen 1400. In the output screen 1400, calculated values ​​of velocity are used as calculation values. The output screen 1400 shows the velocity distribution of the gas-liquid multiphase flow flowing through the flow channel model 1000 at an input time (0.5 s). More specifically, the output screen 1400 shows the velocity distribution of the gas-liquid multiphase flow in a cross section of the flow channel model 1000 cut parallel to the longitudinal direction. The output screen 1400 shows the velocity distribution using shades of color. This allows the velocity distribution of the flow channel model 1000 at any time to be grasped. Note that although the output screen 1400 displays the velocity distribution, it may also display a flow velocity distribution and a concentration distribution.

[0080] 16 is a diagram showing an example of an output screen 1600. In the output screen 1600, the distribution of a gas phase 1007 and a liquid phase 1008 is used as the calculated value of the physical property. The output screen 1600 displays the gas phase 1007 and the liquid phase 1008 in the flow path model 1000 in a distinguishable manner. This allows the distribution of the gas phase 1007 and the liquid phase 1008 in the flow path model 1000 at an arbitrary time step to be grasped.

[0081] The data transfer unit 16 converts the gas-liquid interface information and the second gas-liquid multiphase flow information into a predetermined data format and stores the converted gas-liquid interface information and second gas-liquid multiphase flow information in the memory 40. The predetermined data format is, for example, a data format that can be used in certain software. Alternatively, the predetermined data format may be a list format, or a graph or video data format. The data transfer unit 16 may transfer the gas-liquid interface information and the second gas-liquid multiphase flow information to an external server or another computer.

[0082] The display unit 30 is, for example, a display device such as a liquid crystal display, etc. The display unit 30 displays a display screen of various information output from the output unit 13.

[0083] The memory 40 is, for example, a storage device such as a solid state drive. The memory 40 stores, for example, information output from the data transfer unit 16. The memory 40 stores an information processing program. The memory 40 stores various data necessary for the simulation performed by the information processing system 1.

[0084] 2 is a block diagram showing a detailed configuration of the calculation process unit 12. The calculation process unit 12 includes a previous process acquisition unit 121, a fluid calculation unit 122, a chemical reaction calculation unit 123, a heat transfer calculation unit 124, a diffusion calculation unit 125, a physical property calculation unit 126, a gas-liquid interface calculation unit 127, and an overall calculation unit 128.

[0085] The previous process acquisition unit 121 sets the current time step in the simulation as time step t ｊ Then, at time step t ｊ－１ In the process, the pre-process acquisition unit 121 acquires various calculation values ​​calculated by the fluid calculation unit 122, the chemical reaction calculation unit 123, the heat transfer calculation unit 124, the diffusion calculation unit 125, the physical property calculation unit 126, and the gas-liquid interface calculation unit 127. １ At time step t, the first gas-liquid multiphase flow information, the structure information, and the reaction information acquired by the initial condition acquisition unit 11 are acquired. ｊ In this embodiment, the simulation period is divided into J time steps tｊ It is divided into.

[0086] The fluid calculation unit 122 calculates the time step t １ In the method, the pressure and flow velocity are calculated for each of a plurality of mesh points by applying the initial value of pressure, the initial value of flow velocity, the initial value of viscosity coefficient, the initial value of density, the initial value of gas-liquid interface, structural information, and the initial value of surface tension acquired by the initial condition acquisition unit 11 to the Navier-Stokes equations.

[0087] The fluid calculation unit 122 calculates the time step t ２ Thereafter, the pressure and the flow velocity are calculated for each of a plurality of mesh points by applying the pressure, the flow velocity, the viscosity coefficient, the density, the gas-liquid interface, the structural information, and the surface tension acquired by the previous process acquisition unit 121 to the Navier-Stokes equation. ２ is the time step t １ It is a later time step than

[0088] The chemical reaction calculation unit 123 calculates the time step t １ In the calculation, at least one of the amount of reactant, heat of solution, heat of reaction, and latent heat is calculated based on the initial value of pressure, the initial value of temperature, and at least one of the initial values ​​of compound, chemical reaction formula, reaction rate constant, and molar enthalpy, and the structural information acquired by the initial condition acquisition unit 11. The heat of reaction is calculated when a change occurs in the amount of reactant.

[0089] The chemical reaction calculation unit 123 calculates the time step t ２ Thereafter, based on the pressure, temperature, at least one of the compound, chemical reaction formula, reaction rate constant, and molar enthalpy, and the structural information acquired by the previous process acquisition unit 121, at least one of the amount of reaction component, heat of solution, heat of reaction, and latent heat is calculated for each of multiple mesh points.

[0090] The heat transfer calculation unit 124 calculates the time step t １In this method, at least one of the initial value of the flow rate, the initial value of the amount of reaction components, the initial value of the heat of reaction, and the initial value of the heat of solution, which are acquired by the initial condition acquisition unit 11, the initial value of the constant pressure heat capacity, the initial value of the thermal conductivity, and the structural information are applied to the heat conduction equation to calculate the temperature for each of the multiple mesh points.

[0091] The heat transfer calculation unit 124 calculates the time step t ２ Thereafter, the temperature is calculated for each of the multiple mesh points by applying at least one of the flow rate, amount of reaction components, heat of reaction, and heat of solution acquired by the previous process acquisition unit 121, the constant pressure heat capacity, the thermal conductivity, and the structural information to the heat conduction equation.

[0092] The diffusion calculation unit 125 calculates the time step t １ In the method, at least one of the initial value of the flow rate and the initial value of the amount of reaction components, the initial value of the temperature, the initial value of the constant pressure heat capacity, the initial value of the thermal conductivity, the initial value of the gas-liquid interface and the initial value of the solvation free energy, the initial value of the diffusion coefficient, and the structural information, which are acquired by the initial condition acquisition unit 11, are applied to the diffusion equation to calculate the concentration for each of multiple mesh points.

[0093] The diffusion calculation unit 125 calculates the time step t ２ Thereafter, the concentration is calculated for each of the multiple mesh points by applying to the diffusion equation at least one of the flow rate and the amount of reactive components, the temperature, the constant pressure heat capacity, the thermal conductivity, the gas-liquid interface, and the solvation free energy, the diffusion coefficient, and the structural information acquired by the previous process acquisition unit 121.

[0094] The physical property calculation unit 126 calculates the time step t １ In this calculation, the viscosity coefficient, density, thermal conductivity, constant pressure heat capacity, latent heat of vaporization, diffusion coefficient, and surface tension are calculated for each of multiple mesh points based on the initial value of pressure, initial value of flow velocity, initial value of temperature, initial value of amount of reaction component, initial value of gas-liquid interface, and initial value of concentration acquired by the initial condition acquisition unit 11.

[0095] The physical property calculation unit 126 calculates the time step t ２Thereafter, based on the pressure, flow velocity, temperature, amount of reactive components, gas-liquid interface, and concentration acquired by the fluid calculation unit 122, the viscosity coefficient, density, thermal conductivity, constant pressure heat capacity, latent heat of vaporization, diffusion coefficient, and surface tension are calculated for each of multiple mesh points.

[0096] The gas-liquid interface calculation unit 127 calculates the time step t １ In the method, the gas-liquid interface is calculated by applying the initial value of pressure, the initial value of flow velocity, the initial value of temperature, the initial value of free energy, and / or the initial value of solvation free energy, which are acquired by the initial condition acquisition unit 11, to the Cahn-Hilliard or Allen-Cahn equation.

[0097] The gas-liquid interface calculation unit 127 calculates the time step t ２ Thereafter, the gas-liquid interface is calculated by applying the pressure, flow rate, temperature, free energy, and / or solvation free energy, and structural information acquired by the pre-process acquisition unit 121 to the Cahn-Hilliard or Allen-Cahn equation.

[0098] The overall calculation unit 128 calculates step absolute residuals and step relative residuals of the solutions in each of the fluid calculation unit 122, the chemical reaction calculation unit 123, the heat transfer calculation unit 124, the diffusion calculation unit 125, the physical property calculation unit 126, and the gas-liquid interface calculation unit 127, and determines whether the step absolute residuals and the step relative residuals satisfy the step convergence condition. When the overall calculation unit 128 determines that the step absolute residuals and the step relative residuals of the solutions in each of the fluid calculation unit 122, the chemical reaction calculation unit 123, the heat transfer calculation unit 124, the diffusion calculation unit 125, the physical property calculation unit 126, and the gas-liquid interface calculation unit 127 satisfy the step convergence condition, ｊ at the next time step t ｊ＋１ The step convergence condition is ｊ This is the condition that indicates that the solution has converged.

[0099] On the other hand, when the overall calculation unit 128 determines that the step absolute residual and the step relative residual of the solution in each of the fluid calculation unit 122, the chemical reaction calculation unit 123, the heat transfer calculation unit 124, the diffusion calculation unit 125, the physical property calculation unit 126, and the gas-liquid interface calculation unit 127 do not satisfy the step convergence condition, the overall calculation unit 128 ｊ Continue processing.

[0100] FIG. 3 is a flowchart showing an example of processing performed by the information processing system 1 according to the embodiment.

[0101] (Step S1) The initial condition acquisition unit 11 executes an initial condition acquisition process, the details of which will be described later with reference to FIG.

[0102] (Step S2) The pre-process acquisition unit 121 executes a pre-process acquisition process, the details of which will be described later with reference to FIG.

[0103] (Step S3) The calculation process unit 12 executes a calculation process calculation process, the details of which will be described later with reference to FIGS.

[0104] (Step S4) The overall calculation unit 128 calculates the time step t ｊ The step absolute residual and the step relative residual of the solution at are calculated, and it is determined whether or not the calculated step absolute residual and the step relative residual satisfy the step convergence condition.

[0105] The overall calculation unit 128 calculates the time step t ｊIf the step absolute residue and the step relative residue satisfy step absolute residue < reference absolute residue ε3 and step relative residue < reference relative residue ε4, it may be determined that the step convergence condition is satisfied. Values ​​set in advance by the user are used for the reference absolute residue ε3 and the reference relative residue ε4. Values ​​corresponding to the solutions of the fluid calculation unit 122, the chemical reaction calculation unit 123, the heat transfer calculation unit 124, the diffusion calculation unit 125, the physical property calculation unit 126, and the gas-liquid interface calculation unit 127 are respectively used for the reference absolute residue ε3 and the reference relative residue ε4. The overall calculation unit 128 may determine that the step convergence condition is satisfied if the number of mesh points that do not satisfy step absolute residue < reference absolute residue ε3 is equal to or less than the reference number of mesh points and the number of mesh points that do not satisfy step relative residue < reference relative residue ε4 is equal to or less than the reference number of mesh points. If the overall calculation unit 128 determines that the step convergence condition is not satisfied, the output unit 13 may output failure data indicating the mesh points at which the calculation failed. For example, the output unit 13 may display the mesh points for which the calculation has failed on the display screen of the flow channel model 1000 .

[0106] If the step convergence condition is satisfied (YES in step S4), the process proceeds to step S7, and if the step convergence condition is not satisfied (NO in step S4), the process proceeds to step S5.

[0107] (Step S5) The overall calculation unit 128 determines whether or not the processing of steps S3 and S4 has been performed K times. If the processing has not been performed K times (NO in step S5), the process proceeds to step S6. If the processing has been performed K times (YES in step S5), the process proceeds to step S8. If further calculations are performed for K, the time step t ｊ A predetermined value is adopted that is considered to be unlikely to result in a solution at .times. ...

[0108] (Step S6) The overall calculation unit 128 changes the convergence parameter and returns the process to step S3. The convergence parameter will be described later. In this case, the time step t ｊ The calculations at are performed.

[0109] (Step S7) The output unit 13 outputs the time step tｊ The simulation results in the calculation unit 1200 are output to the display unit 30 or the like. The simulation results include the pressure and flow velocity calculated by the fluid calculation unit 122. The simulation results include at least one of the amounts of reactant components, heat of solution, heat of reaction, and latent heat calculated by the chemical reaction calculation unit 123. The simulation results include the temperature calculated by the heat transfer calculation unit 124. The simulation results include the concentration calculated by the diffusion calculation unit 125. The simulation results include the viscosity coefficient, density, thermal conductivity, constant pressure heat capacity, latent heat of vaporization, diffusion coefficient, and surface tension calculated by the physical property calculation unit 126. The simulation results include the gas-liquid interface calculated by the gas-liquid interface calculation unit 127. The simulation results include the output screen 1300, the output screen 1400, and the output screen 1600.

[0110] (Step S8) The overall calculation unit 128 calculates the time step t ｊ is the final time step t Ｊ Determine whether the time step t ｊ is the final time step t Ｊ If the time step t ｊ is the final time step t Ｊ If it has not reached (NO in step S8), the process proceeds to step S9.

[0111] (Step S9) The overall calculation unit 128 calculates the time step t ｊ at the next time step t ｊ＋１ Then, the process returns to step S2. ｊ＋１ A simulation is performed at

[0112] FIG. 4 is a flowchart illustrating an example of the initial condition acquisition process.

[0113] (Step S101) The initial condition acquisition unit 11 acquires first gas-liquid multiphase flow information. The first gas-liquid multiphase flow information includes an initial value of pressure, an initial value of flow velocity, an initial value of temperature, an initial value of concentration, an initial value of constant pressure heat capacity, an initial value of thermal conductivity, an initial value of viscosity coefficient, an initial value of density, an initial value of gas-liquid interface, an initial value of diffusion coefficient, a type of substance, an initial value of free energy, an initial value of surface tension, and an initial value of solvation free energy.

[0114] (Step S102) The initial condition acquisition unit 11 acquires structural information, which includes the flow channel model 1000 and boundary conditions.

[0115] (Step S103) The initial condition acquiring unit 11 acquires reaction information, which includes a chemical reaction formula, a reaction rate constant, molar enthalpy, initial values ​​of latent heat of vaporization, and initial values ​​of heat of reaction.

[0116] FIG. 5 is a flowchart showing an example of the pre-processing.

[0117] (Step S200) The previous process acquisition unit 121 acquires the time step t ｊ is the time step t １ At time step t ｊ is the time step t １ If the time step t ｊ is the time step t １ If not (NO in step S200), that is, if the time step t ｊ is the time step t ２ If so, the process proceeds to step S201.

[0118] (Step S201) The previous process acquisition unit 121 acquires the time step t ｊ－１ The pressure and flow velocity calculated in the above are obtained.

[0119] (Step S202) The previous process acquisition unit 121 calculates the time step t ｊ－１ The amounts of reaction components, heat of dissolution, heat of reaction, and latent heat calculated in the above are obtained.

[0120] (Step S203) The pre-process acquisition unit 121 calculates the time step t ｊ－１ The calculated temperature is obtained.

[0121] (Step S204) The pre-process acquisition unit 121 calculates the time step t ｊ－１ The calculated concentration is obtained.

[0122] (Step S205) The pre-process acquisition unit 121 calculates the time step t ｊ－１ The viscosity coefficient, density, thermal conductivity, constant pressure heat capacity, latent heat of vaporization, diffusion coefficient, and surface tension calculated in the above are acquired.

[0123] (Step S206) The previous process acquisition unit 121 calculates the time step t ｊ－１ The gas-liquid interface calculated by is obtained.

[0124] (Step S207) The previous process acquisition unit 121 acquires the first gas-liquid multiphase flow information, the structure information, and the reaction information from the initial condition acquisition unit 11.

[0125] Fig. 6 is a flowchart showing the details of the calculation step calculation process. Fig. 7 is a flowchart showing the details of the calculation step calculation process continued from Fig. 6. Fig. 8 is a flowchart showing the details of the calculation step calculation process continued from Fig. 7.

[0126] (Step S301) The fluid calculation unit 122 applies the pressure p, flow velocity u, viscosity coefficient, density, gas-liquid interface, structural information, and surface tension acquired by the previous process acquisition unit 121 to the Navier-Stokes equation and the continuity equation to calculate the pressure p and flow velocity u at the next time for each mesh point.

[0127] The finite element method or boundary element method is used to calculate the behavior of a flow path using CAE (computer-aided engineering). In the finite volume method or finite difference method, the calculation domain is represented by a combination of polygons such as triangles, quadrilaterals, or pentagons. Such polygons are called meshes. The vertices of the mesh are called mesh points. The calculation process unit 12 basically calculates each equation at each mesh point. For sub-mesh points between mesh points, the calculation process unit 12 can interpolate using the calculated values ​​of surrounding mesh points.

[0128] Assuming that the fluid is a Newtonian fluid, the Navier-Stokes equations are expressed by the following equation (1), and the continuity equation is expressed by equation (2).

[0129]

[0130] μ is the shear viscosity, χ is the bulk viscosity, ρ is the density, σ is the surface tension, and ε ２ ｐｆ is a parameter included in the phase field (PF) method described later, and φ is a PF variable. The fluid calculation unit 122 calculates the pressure p and flow velocity u at the next time by solving the Navier-Stokes equations using the flow path model 1000 and boundary conditions included in the structural information as constraints.

[0131] The fluid calculation unit 122 may solve the Navier-Stokes equations using various methods, such as the Reynolds Averaged Model (RANS), Large Eddy Simulation (LES), or Direct Numerical Simulation (DNS).

[0132] When using RANS, the fluid calculation unit 122 does not solve the Navier-Stokes equations directly, but models the Navier-Stokes equations by making certain assumptions, and then solves the Navier-Stokes equations using this model. As models, yPlus, L-VEL, k-ε, reanalyzable k-ε, reanalyzable k-ω, SST, Spaart-Allmaras, and v2-f can be adopted. Furthermore, SST, Spaart-Allmaras, or v2 may be adopted as models of the Navier-Stokes equations, in which the turbulence parameters required for RANS calculations may vary depending on the distance from the wall, etc.

[0133] The fluid calculation unit 122 may obtain a solution using a RANS model, which solves the Navier-Stokes equations under certain assumptions for both mesh points and sub-mesh points, or may obtain a solution using an LES model or direct numerical simulation (DNS), which has a higher computational load than a RANS model. For example, the LES model may be the Smagnorinsky model, the dynamic Smagnorinsky model, the RBVM model, the RBMVWV model, or a scale similarity model. More preferably, the LES model may be the dynamic Smagnorinsky model or the RBVMWV model, which do not use experimental results.

[0134] In general, the appearance of the Navier-Stokes equations can change depending on several assumptions, such as the assumption of incompressibility (density does not change), the assumption of compressibility (density changes), and the assumption that the viscosity coefficient is constant or changes. The fluid calculation unit 122 may adopt any of these assumptions. For example, the fluid calculation unit 122 may ignore the bulk modulus χ, which is relatively small, and adopt the assumption that only the shear viscosity coefficient μ depends on the temperature T and pressure p, or the assumption that only the shear viscosity coefficient μ depends on the gas-liquid concentration.

[0135] (ρ-ρ) in formula (1) ｌｉｑ ) g / ρ is the buoyancy term of the air bubble, ε ２ ｐｆ φ∇∇ ２ φ−(2√2 / (3σε ｐｆ ))φ(3φ ２+1) ∇φ is a term related to surface tension. The surface tension σ is calculated based on the temperature or dissolution calculated based on reaction information as described below. The surface tension of the entire gas phase can be calculated by calculating the term related to surface tension using the value of surface tension σ. The gas-liquid interface can be calculated so that the calculated surface tension of the entire gas phase matches the external pressure p. In particular, the effect of the surface tension of the entire gas phase is important in determining the gas-liquid interface. ρ ｌｉｑ is the density ρ of the liquid. The fluid calculation unit 122 uses the values ​​acquired by the initial condition acquisition unit 11 in the initial stage of the simulation, and after repeated calculations, calculates the pressure p and the flow velocity u using the values ​​acquired by the previous process acquisition unit 121. Note that the surface tension σ may be calculated based on the temperature and dissolution calculated based on the reaction information.

[0136] The boundary conditions are constituted by the inlet 1005, the outlet 1006, the wall surface of the flow path 1002, and the wall surface of the heat source 1003. The fluid calculation unit 122 mainly determines the flow velocity u for the inlet 1005. On the other hand, for the outlet 1006, the fluid calculation unit 122 may set the pressure p at the time of fluid outflow to p0, and if the pressure p is smaller than p0, set the pressure p to p0 to prevent backflow. The fluid calculation unit 122 may perform calculations assuming that there is no slippage of the gas-liquid multiphase flow for the wall surface of the flow path 1002. In this case, the fluid calculation unit 122 may set the flow velocity u at the wall surface to u=0.

[0137] (Step S302) The fluid calculation unit 122 determines whether the solution of step S301 satisfies the convergence condition. The fluid calculation unit 122 calculates the absolute residual and the relative residual of the nth solution for each mesh point. n is an index that specifies the number of calculations of step S301. The absolute residual is specified by |left-hand side value - right-hand side value|. n indicates the number of times the equation is solved, that is, the number of times the processing of step S301 is executed. The solution here is a solution to the Navier-Stokes equation, that is, either the pressure p or the flow velocity u. The left-hand side value refers to the value of the left-hand side when the nth solution is substituted into the equation. The right-hand side value refers to the value of the right-hand side when the nth solution is substituted into the equation. The relative residual is specified by "absolute residual / maximum value of absolute residual at all mesh points."

[0138] The fluid calculation unit 122 may determine that the convergence condition is satisfied if the absolute residue is smaller than the reference absolute residue ε1 and the relative residue is smaller than the reference relative residue ε2 at all mesh points. The reference absolute residue ε1 and the reference relative residue ε2 are values ​​preset by the user. The fluid calculation unit 122 may also determine that the convergence condition is satisfied if the number of mesh points that do not satisfy the condition absolute residue < reference absolute residue ε1 is equal to or less than a reference number and the number of mesh points that do not satisfy the condition relative residue < reference relative residue ε2 is equal to or less than a reference number. If the overall calculation unit 128 determines that the convergence condition is not satisfied, the output unit 13 outputs failure data indicating the mesh points at which the calculation failed. For example, the output unit 13 may display the mesh points at which the calculation failed on the display screen of the flow path model 1000. The same applies to steps S306, S310, S314, S318, and S322.

[0139] If the convergence condition is satisfied (YES in step S302), the process proceeds to step S305; if the convergence condition is not satisfied (NO in step S302), the process proceeds to step S303.

[0140] (Step S303) The fluid calculation unit 122 determines whether the number of calculations in step S301 has reached N. If the number of calculations has reached N (YES in step S303), the process proceeds to step S7 in FIG. 3. In this case, the simulation results up to the middle of the simulation period are output. If the number of calculations has not reached N (NO in step S303), the process proceeds to step S304. N is a predetermined value.

[0141] (Step S304) The fluid calculation unit 122 changes the convergence parameter, and the process returns to step S301. The convergence parameter will be described later.

[0142] (Step S305) The chemical reaction calculation unit 123 calculates the amounts of reaction components, the heat of solution R based on the pressure p, the temperature T, the compound, the chemical reaction formula, at least one of the reaction rate constant and the molar enthalpy, and the structural information acquired by the previous process acquisition unit 121. ｄｉｖ , reaction heat R ｒｅａｃ , and latent heat R ｐｃis calculated for each of the mesh points.

[0143] The number of compounds constituting each of the liquid phase and the gas phase may be any number of types as long as it is one or more. It is assumed that compounds in each phase can be mixed if a sufficient time is allowed. However, if a part of the liquid phase or the gas phase is not mixed and a new phase is visually recognized, the chemical reaction calculation unit 123 changes the type of the PF variable φ as described later, for example, φ Ａ , φ Ｂ In this way, we can increase the number of elements and calculate the solution.

[0144] A chemical reaction is composed of one or more chemical reaction formulas. A chemical reaction formula may include equilibrium in part or all of the chemical reaction formulas. A chemical reaction formula may involve a chemical reaction in at least one of the gas phase and the liquid phase. For the sake of convenience, the following description will be given assuming that only the gas phase undergoes a chemical reaction.

[0145] However, this is just an example, and the liquid phase may be vaporized for a certain component k in the liquid phase. ｋ ｔｏｔ is the total molar amount of component k, n ｋ ｌｉｑ is the liquid phase molar amount of component k, n ｋ ｇａｓ When is the molar amount of component k when vaporized, the chemical reaction formula is expressed by formula (3).

[0146]

[0147] In this disclosure, chemical reactions include not only chemical reactions such as decomposition and polymerization, but also dissolution reactions of gas phase into liquid phase. However, in calculations, it is easier to understand if the two are treated separately. Therefore, this disclosure is referred to as the heat of solution R ｄｉｖ and reaction heat R ｒｅａｃ Generally, the calorific value R ｔｏｔ is the reaction heat R ｒｅａｃ , the heat of solution is R ｄｉｖ , latent heat R ｐｃ Then, it can be written as in equation (4).

[0148]

[0149] In the present disclosure, chemical reactions include dissolution reactions, and therefore the heat of reaction Rｒｅａｃ = 0 and the heat of solution R ｄｉｖ = 0 cannot be true at the same time.

[0150] Latent heat R ｐｃ is the latent heat of vaporization of component k at temperature T. ｐｃ Then, it is expressed by the formulas (5) and (6).

[0151]

[0152] If the free energy F shown in the equation (31) etc. described later is used, the latent heat R ｐｃ can be easily determined. Note that a function satisfying ζ(φ) is selected such that ζ(-1)=0, ζ(0)=0.5, ζ(1)=1, -dζ(-1) / dφ=0, and -dζ(1) / dφ=0. Therefore, the chemical reaction calculation unit 123 calculates the latent heat R ｐｃ Calculate.

[0153] The physical property values, chemical potential, and equation of state are obtained by experimental values, quantum chemical calculations, and microsimulations such as molecular dynamics, etc. Alternatively, the physical property values, chemical potential, and equation of state may be obtained by model calculations based on certain assumptions.

[0154] When a chemical reaction is taken into consideration, the amount of reactants and the heat of reaction R ｒｅａｃ In the s-th chemical reaction equation of the gas phase component, the coefficient of component k is ν ｓｋ Then, for example, it is expressed by equation (7).

[0155]

[0156] The rate constant of the forward reaction is k ｆ ＡＢ , the rate constant of the reverse reaction is k ｒ ＡＢ In this equation, for example, the rates of the reaction components are expressed by equations (8), (9), (10), and (11).

[0157]

[0158] In the above formula, the amount of reaction components refers to the amounts [C] and [D] of components C and D that change due to the chemical reaction of components A and B, and the amounts [A] and [B] of components A and B that change due to the chemical reaction. These component amounts can be expressed in various ways, such as molar amount, volume amount, and mass, but the molar concentration, which is the molar amount divided by the volume, is easy to handle.

[0159] The chemical reaction calculation unit 123 may calculate the amounts of reaction components of the compound using formulas (8) to (11). Components A to D are examples of compounds or substances. Generally, a chemical reaction is composed of multiple reaction formulas, and therefore the above formulas can occur for as many amounts of components.

[0160] Reaction rate constant k ｆ ｉ and reaction rate constant k ｒ ｉ Assume that the temperature follows the Arrhenius equation. In this case, the temperature is T, the gas constant is R, and the activation energies of the forward and reverse reactions determined in advance through experiments or simulations are E. ｆ , E ｒ , the frequency factors of forward and reverse reactions are A ｆ , A ｒ , the temperature exponents of the forward and reverse reactions are n ｆ , n ｒ Then, the reaction rate constant k ｆ ｉ and reaction rate constant k ｒ ｉ are expressed by equations (12) and (13).

[0161]

[0162] This gives the reaction rate constant k ｆ ｉ and reaction rate constant k ｒ ｉ can be determined by the temperature T. If the values ​​of the activation energy, frequency factor, and temperature exponent are known in advance, those values ​​are used. If the reaction rate constant is known directly, that value is used. These values ​​are acquired by the initial condition acquisition unit 11.

[0163] The molar enthalpy hr generated by the reaction is expressed by equations (14) and (15).

[0164]

[0165] React is the reactant (A and B in formula (7)), and prod is the compound produced (C and D in formula (7)). ｋ is the molar enthalpy of each component k in the forward and reverse reactions in the i-th chemical reaction equation. ｌ is the molar enthalpy of each component l in the forward and reverse reactions in the i-th chemical reaction equation. ｆ ｉ is the sum of the molar enthalpies of the forward reactions. ｒ ｉ is the sum of the molar enthalpies of the reverse reactions. The values ​​of react and prod are acquired by the initial condition acquisition unit 11. The sum of the molar enthalpies of the forward reactions, h ｆ ｉ , the sum of the molar enthalpies of the reverse reactions h ｒ ｉ is acquired by the initial condition acquisition unit 11, the chemical reaction calculation unit 123 calculates the sum of the molar enthalpies of the forward reactions h without using equations (14) and (15). ｆ ｉ , the sum of the molar enthalpies of the reverse reactions h ｒ ｉ h ｋ , h ｌ is a coefficient a of a function shown in the following formula (X), such as a NASA polynomial, in the initial condition acquisition unit 11. １，ｍ ~a ６，ｍ where T is the temperature and R is the gas constant. m is substituted with k or l. prod is given by v ｓｋ ・h ｋ is applicable.

[0166]

[0167] By using these equations, the reaction heat R ｒｅａｃ is written as follows:

[0168]

[0169] Assuming that component k in the gas phase dissolves in the liquid phase, the heat of dissolution of component k is Δh ｄｉｖ Then, the heat of solution R ｄｉｖis the latent heat R shown in equations (5) and (6). ｐｃ Similarly, it is expressed by the formula (17). Therefore, the chemical reaction calculation unit 123 calculates the heat of dissolution R using the formula (17) from the free energy F shown in the formula (31) etc. ｄｉｖ Calculate.

[0170]

[0171] From the above, the reaction heat R ｒｅｃｔ , latent heat R ｐｃ , and the heat of solution R ｄｉｖ The amount of heat generated by R ｔｏｔ is calculated.

[0172] The chemical reaction calculation unit 123 calculates the amounts of reaction components and the reaction heat R ｒｅａｃ and the heat of solution R ｄｉｖ In this case, the chemical reaction calculation unit 123 may assume that no chemical reaction occurs on the wall surface of the flow channel 1002. That is, the chemical reaction calculation unit 123 calculates the reaction heat R ｒｅｃｔ and heat of solution R ｄｉｖ R ｒｅａｃ = 0, R ｄｉｖ = 0.

[0173] As described above, the chemical reaction calculation unit 123 calculates, as an example, the reaction rate constant k ｆ ＡＢ , k ｒ ＡＢ The amounts of reaction components are calculated by substituting the above and calculating the amounts of components [A], [B], [C], and [D].

[0174] The compounds correspond to components A to D shown in formulas (7) to (11). The chemical reaction formula corresponds to, for example, formula (7). The reaction rate constant is, for example, the reaction rate constant k ｆ ＡＢ , k ｒ ＡＢ Here, the chemical reaction calculation unit 123 calculates the reaction rate constant k ｆ ＡＢ , k ｒ ＡＢ If is unknown, the reaction rate constant k can be obtained by substituting the temperature T into equations (12) and (13).ｆ ＡＢ , k ｒ ＡＢ can be calculated.

[0175] The chemical reaction calculation unit 123 calculates the sum of the molar enthalpies of the forward reaction h ｆ ｉ , the sum of the molar enthalpies of the reverse reactions h ｒ ｉ By substituting, the reaction heat R ｒｅａｃ Calculate.

[0176] The molar enthalpy is the sum of the molar enthalpies of the forward reactions in equations (14) and (15), h ｆ ｉ , the sum of the molar enthalpies of the reverse reactions h ｒ ｉ is applicable.

[0177] The total molar enthalpy of the forward reaction is h ｆ ｉ and the sum of the molar enthalpies of the reverse reactions h ｒ ｉ When the molar enthalpies h and h are unknown, the chemical reaction calculation unit 123 calculates the sum of the molar enthalpies h of the forward reactions using equations (14) and (15). ｆ ｉ , the sum of the molar enthalpies of the reverse reactions h ｒ ｉ can be calculated.

[0178] The chemical reaction calculation unit 123 calculates the latent heat R by substituting the temperature T and the free energy F into the formula (5). ｐｃ The free energy F is expressed by the equation (31) described below. The pressure p is used when determining the free energy F.

[0179] The chemical reaction calculation unit 123 calculates the heat of solution R by substituting the temperature T and the free energy F into equation (17). ｄｉｖ The free energy F is expressed by the equation (31) described below. The pressure p is used when determining the free energy F.

[0180] (Step S306) The chemical reaction calculation unit 123 determines whether the calculation result of step S305 satisfies the convergence condition. As in step S302, the chemical reaction calculation unit 123 calculates the absolute residue and the relative residue of the n-th solution for each mesh point. The solution here is, for example, the amount of reaction component calculated using the equations (Equations (8) to (11)).

[0181] The chemical reaction calculation unit 123 may determine that the convergence condition is satisfied if the absolute residue < the reference absolute residue ε1 and the relative residue < the reference relative residue ε2. Values ​​set in advance by the user are adopted as the reference absolute residue ε1 and the reference relative residue ε2. Values ​​suitable for the solution of the chemical reaction calculation unit 123 are adopted as the reference absolute residue ε1 and the reference relative residue ε2. The chemical reaction calculation unit 123 may determine that the convergence condition is satisfied if the number of mesh points that do not satisfy the condition absolute residue < the reference absolute residue ε1 is equal to or less than a reference number and the number of mesh points that do not satisfy the condition relative residue < the reference relative residue ε2 is equal to or less than a reference number.

[0182] If the convergence condition is satisfied (YES in step S306), the process proceeds to step S309; ​​if the convergence condition is not satisfied (NO in step S306), the process proceeds to step S307.

[0183] The processes in steps S307 and S308 are the same as those in steps S303 and S304.

[0184] (Step S309) The heat transfer calculation unit 124 calculates the flow rate u, the amount of reaction components, and the heat of dissolution R ｄｉｖ , and the reaction heat R ｒｅａｃ At least one of the above and constant pressure heat capacity C ｐ and thermal conductivity k Ｔ and the structural information are applied to the heat conduction equation to calculate the temperature T for each of the plurality of mesh points.

[0185] The heat transfer calculation unit 124 calculates the latent heat R ｐｃ and the heat of solution R ｄｉｖ and reaction heat R ｒｅａｃ The heat generation amount R is the sum of ｔｏｔ and the flow velocity u calculated by the fluid calculation unit 122 may be substituted into the heat conduction equation shown in equation (18) below.

[0186] V ｌｉｑ is the volume ratio of the liquid phase in the fluid, V ｇａｓ Let be the volume ratio of the gas phase in the fluid. Also, let k be a component in the gas phase, and l be a component in the liquid phase. Then, let C be the constant pressure heat capacity, thermal conductivity, density, and molar ratio of component k, respectively. ｐ＿ｇａｓ ｋ , k Ｔ＿ｇａｓ ｋ , ρ ｇａｓ ｋ , m ｇａｓ ｋ The constant pressure heat capacity, thermal conductivity, density, and molar ratio of component l are C ｐ＿ｌｉｑ ｌ , k Ｔ＿ｌｉｑ ｌ , ρ ｌｉｑ ｌ , m ｌｉｑ ｌ In this case, the heat conduction equation is written as equation (18).

[0187]

[0188] C ｐ , k Ｔ , ρ are expressed by equations (18) to (21).

[0189]

[0190] Volume ratio V ｌｉｑ , V ｇａｓ is expressed by Equation (22). In a system to which the phase field method is applied, the gas-liquid interface is assumed to be smoothly connected and symmetric about the gas-liquid interface. Therefore, φ is written in the form of tanh as shown in Equation (23).

[0191] m ｇａｓ ｋ corresponds to the amount of reactant. For example, if components A and C are gases, m ｇａｓ ａ = [A], and m ｇａｓ ｃ = [C].

[0192]

[0193] More specifically, the heat transfer calculation unit 124 calculates V using the gas-liquid interface acquired by the previous process acquisition unit 121. ｇａｓ, V ｌｉｑ The heat transfer calculation unit 124 determines the concentration c of the gas component i. ｇａｓ ｋ , and the concentration of component j in the liquid, c ｌｉｑ ｊ are given as initial values, it is possible to track the change in temperature T over time using equation (18).

[0194] Molar ratio m of component k in the gas phase ｇａｓ ｋ is the molar concentration of component k in the gas phase, ｇａｓ ｋ Then, the molar concentration m of component l in the liquid phase is expressed by equation (24). ｌｉｑ ｌ is the molar concentration of component k dissolved in the liquid, ｇａｓ＿ｄｉｖ ｋ Then, it is expressed by equation (25).

[0195] For example, if components A and C are gases and component A dissolves in components B and D, then ｇａｓ＿ｄｉｖ ｋ =C ｇａｓ＿ｄｉｖ Ａ Then, Σ ｌ C ｌｉｑ ｌ =C ｇａｓ Ｃ +C ｇａｓ Ａ Then, C ｇａｓ Ｃ = [C], and C ｇａｓ Ａ +C ｇａｓｄｉｖ Ａ = [A].

[0196]

[0197] In practice, the heat transfer calculation unit 124 calculates the time step t ｊ In this case, φ calculated by the gas-liquid interface calculation unit 127 and the concentration c calculated by the diffusion calculation unit 125 are ｇａｓ＿ｄｉｖ ｋ However, in this case, the heat transfer calculation unit 124 calculates the temperature T by simultaneously calculating the time step t ｊ－１ The temperature T is calculated using the calculation results in (a).

[0198] The boundary conditions of the heat conduction equation mainly include the inflow of fluid into the inlet 1005, the outflow of fluid from the outlet 1006, and heat conduction from the wall surface of the flow path 1002. The inflow of fluid into the inlet 1005 is mainly determined by the heat quantity Q or the temperature T. The outflow of fluid from the outlet 1006 is determined as un*Q=0, where un is the flow velocity u perpendicular to the wall surface. Heat conduction from the wall surface is determined as un*Q=0 in the case of adiabatic conditions, and is determined by fixing the temperature T in the case of a constant temperature.

[0199] (Step S310) The heat transfer calculation unit 124 determines whether the solution of step S309 satisfies the convergence condition. As in step S302, the heat transfer calculation unit 124 calculates the absolute residual and the relative residual of the n-th solution for each mesh point. The solution here is the solution of the heat conduction equation shown in formula (18), i.e., the temperature T.

[0200] The heat transfer calculation unit 124 may determine that the convergence condition is satisfied if the absolute residue is less than the reference absolute residue ε1 and the relative residue is less than the reference relative residue ε2. Values ​​set in advance by the user are adopted as the reference absolute residue ε1 and the reference relative residue ε2. Values ​​suitable for the solution obtained by the heat transfer calculation unit 124 are adopted as the reference absolute residue ε1 and the reference relative residue ε2. The heat transfer calculation unit 124 may determine that the convergence condition is satisfied if the number of mesh points that do not satisfy the condition absolute residue < reference absolute residue ε1 is equal to or less than a reference number and the number of mesh points that do not satisfy the condition relative residue < reference relative residue ε2 is equal to or less than a reference number.

[0201] If the convergence condition is satisfied (YES in step S310), the process proceeds to step S313; if the convergence condition is not satisfied (NO in step S310), the process proceeds to step S311.

[0202] The processes in steps S311 and S312 are the same as those in steps S303 and S304 in FIG.

[0203] (Step S313) The diffusion calculation unit 125 calculates at least one of the flow velocity u and the amount of reactant, the temperature T, and the constant pressure heat capacity C acquired by the previous process acquisition unit 121. ｐ and thermal conductivity k ＴThe concentration is calculated for each of a plurality of mesh points by applying at least one of the gas-liquid interface and the solvation free energy, the diffusion coefficient D, and the structural information (the flow path model 1000 and the boundary conditions) to the diffusion equation.

[0204] When a component k in the gas phase dissolves in the liquid phase, the concentration c ｇａｓ＿ｄｉｖ ｋ The time evolution equation for is expressed by equation (26). The diffusion calculation unit 125 calculates the free energy F ｋ , mobility M φ , density ρ ｇａｓ ｋ , and the diffusion coefficient D ｇａｓ＿ｄｉｖ ｋ , D ｇａｓ＿ｄｉｖ＿Ｔ ｋ is substituted into the diffusion equation shown in equation (26), the concentration c ｇａｓ＿ｄｉｖ ｋ , c ｇａｓ ｋ Calculate the following: ｇａｓ ｋ is the concentration of component k in the gas phase. The gas-liquid interface is the PF variable φ shown in equation (26). Note that the structural information is used as a constraint when solving equation (26).

[0205]

[0206] The diffusion calculation unit 125 calculates the concentration c ｇａｓ ｋ , c ｇａｓ＿ｄｉｖ ｋ For example, when component A dissolves, [A] = C ｇａｓ ｋ +C ｇａｓ＿ｄｉｖ ｋ In reality, the heat of solution R occurs when a component i in the gas phase dissolves in the liquid phase. ｄｉｖ In some cases, the change in ∇T due to the ∇T does not change significantly. In this case, the diffusion calculation unit 125 can ignore the ∇T term.

[0207] In the phase-field (PF) method, an equilibrium state is defined as a state in which the chemical potentials of the gas phase and the liquid phase are equal. In other words, when the gas phase is a solute, the chemical potential of the gas phase μ ｒｅｆ (gas) and the liquid phase chemical potential μ ｒｅｆ Regarding (liq), equation (27) holds.

[0208]

[0209] However, locally, the chemical potentials of the liquid phase and the gas phase differ depending on the dissolution of the liquid and gas phases. Therefore, the diffusion calculation unit 125 uses the difference between the liquid phase chemical potential and the gas phase chemical potential as a driving force, and considers concentration diffusion in the second term on the right side of equation (26) (the following equation) so as to reduce the difference.

[0210]

[0211] (Step S314) The diffusion calculation unit 125 determines whether the calculation result of step S313 satisfies the convergence condition. As in step S302, the diffusion calculation unit 125 calculates the absolute residual and the relative residual of the n-th solution for each mesh point. The solution here is the solution of the diffusion equation shown in equation (26), i.e., the concentration.

[0212] The diffusion calculation unit 125 may determine that the convergence condition is satisfied if the absolute residue is smaller than the reference absolute residue ε1 and the relative residue is smaller than the reference relative residue ε2. Values ​​set in advance by the user are adopted as the reference absolute residue ε1 and the reference relative residue ε2. Values ​​suitable for the solution obtained by the diffusion calculation unit 125 are adopted as the reference absolute residue ε1 and the reference relative residue ε2. The diffusion calculation unit 125 may determine that the convergence condition is satisfied if the number of mesh points that do not satisfy the condition "absolute residue < reference absolute residue ε1" is equal to or less than a reference number and the number of mesh points that do not satisfy the condition "relative residue < reference relative residue ε2" is equal to or less than a reference number.

[0213] If the convergence condition is satisfied (YES in step S314), the process proceeds to step S317; if the convergence condition is not satisfied (NO in step S314), the process proceeds to step S315.

[0214] The processes in steps S315 and S316 are the same as those in steps S303 and S304 in FIG.

[0215] (Step S317) The physical property calculation unit 126 calculates the viscosity coefficient, density ρ, and thermal conductivity k based on the pressure p, flow velocity u, temperature T, amount of reaction component, gas-liquid interface, and concentration. Ｔ and constant pressure heat capacity C ｐ and the latent heat of vaporization Δh ｐｃ The diffusion coefficient D and the surface tension σ are calculated for each of the mesh points.

[0216] The pressure p and flow velocity u are calculated by the fluid calculation unit 122. The temperature T is calculated by the heat transfer calculation unit 124. The amounts of reactants are calculated by the chemical reaction calculation unit 123 (for example, [A] to [D]). The gas-liquid interface is calculated by the PF variable φ calculated by the gas-liquid interface calculation unit 127. The concentration is calculated by the diffusion calculation unit 125. ｇａｓ＿ｄｉｖ ｋ , c ｇａｓ ｋ is used.

[0217] For example, the physical property calculation unit 126 calculates the viscosity coefficient, density ρ, thermal conductivity k according to the temperature T and pressure p. Ｔ , constant pressure heat capacity C ｐ , latent heat of vaporization Δh ｐｃ , the diffusion coefficient D, and the surface tension σ are associated with each other, and the viscosity coefficient, density ρ, thermal conductivity k are determined using a property list in which the values ​​are determined in advance by experiment. Ｔ , constant pressure heat capacity C ｐ , latent heat of vaporization Δh ｐｃ , and the diffusion coefficient D may be determined.

[0218] Fig. 18 is a diagram showing an example of a physical property list. The physical property list shown in Fig. 18 shows the physical property list for water (all liquid) at a temperature T of 283 K and a pressure p of 0.1 to 5.0 MPa. This physical property list shows the density ρ, molar enthalpy, constant pressure heat capacity C according to the given temperature T and pressure p. ｐ , thermal conductivity k Ｔ , the viscosity coefficient, and the diffusion coefficient D are determined.

[0219] 19 is a diagram showing another example of a property list. The property list shown in FIG. 19 shows the property list of water (gas and liquid coexist) at temperatures T of 283 to 363 K. In this property list, the vapor pressure and latent heat of vaporization Δh corresponding to a given temperature T are ｐｃ The physical property calculation unit 126 determines the viscosity coefficient, density ρ, thermal conductivity k using the physical property lists shown in FIGS. Ｔ , constant pressure heat capacity C ｐ , latent heat of vaporization Δh ｐｃ , and the diffusion coefficient D may be determined.

[0220] The physical property calculation unit 126 uses some kind of model to calculate the viscosity coefficient, density ρ, thermal conductivity k Ｔ , constant pressure heat capacity C ｐ , latent heat of vaporization Δh ｐｃ , the diffusion coefficient D, and the surface tension σ. The physical property calculation unit 126 may change the temperature T and the pressure p by using quantum chemistry, molecular dynamics, etc., and determine the viscosity coefficient, density ρ, thermal conductivity k Ｔ , constant pressure heat capacity C ｐ , latent heat of vaporization Δh ｐｃ , the diffusion coefficient D, and the surface tension σ. Ｔ , constant pressure heat capacity C ｐ , latent heat of vaporization Δh ｐｃ , diffusion coefficient D, and surface tension σ are fitted with functions using pressure p, flow velocity u, temperature T, gas-liquid interface (PF variable φ), concentration, etc. to obtain the final viscosity coefficient, density ρ, thermal conductivity k Ｔ , constant pressure heat capacity C ｐ , latent heat of vaporization Δh ｐｃ , diffusion coefficient D, and surface tension σ may be calculated. In particular, it has been found that surface tension has a great influence on identifying the gas-liquid interface. Surface tension is easily affected by temperature and dissolution calculated based on reaction information. Furthermore, as described above, the surface tension σ obtained individually is used to calculate the surface tension of the entire gas phase, and the gas-liquid interface can be calculated so that the calculated surface tension of the entire gas phase matches the external pressure, so the surface tension σ is important. Furthermore, for example, the viscosity coefficient of a certain component i in the gas phase can be calculated as μ ｇａｓ ｉ Let the viscosity coefficient of a certain component j in the liquid phase be μ ｌｉｑｊ The fitting coefficient for the temperature T of the gas phase component i is set as ω ｈｇａｓ The fitting coefficient for the density ρ of the gas phase component i is set as ω ｗ＿ｌｉｑ For the liquid phase component j, the fitting coefficient for the temperature T is ω ｗ＿ｌｉｑ For the liquid phase component j, the fitting coefficient for the density ρ is ξ ｒ＿ｌｉｑ In this case, the physical property calculation unit 126 may calculate the viscosity coefficient μ′ using equation (28).

[0221]

[0222] (Step S318) The physical property calculation unit 126 determines whether the calculation result of step S317 satisfies the convergence condition. As in step S302, the physical property calculation unit 126 calculates the absolute residual and the relative residual of the nth solution for each mesh point. The solution here is one or more of the viscosity coefficient, density, thermal conductivity, heat capacity at constant pressure, latent heat of vaporization, and diffusion coefficient. Note that if the physical property calculation unit 126 determines the viscosity coefficient, density, thermal conductivity, heat capacity at constant pressure, latent heat of vaporization, and diffusion coefficient using a physical property list, equations are not used, and the determination process of step S318 is omitted.

[0223] The physical property calculation unit 126 may determine that the convergence condition is satisfied if the absolute residue is smaller than the reference absolute residue ε1 and the relative residue is smaller than the reference relative residue ε2. Values ​​set in advance by the user are adopted as the reference absolute residue ε1 and the reference relative residue ε2. Values ​​suitable for the solution of the physical property calculation unit 126 are adopted as the reference absolute residue ε1 and the reference relative residue ε2. The physical property calculation unit 126 may determine that the convergence condition is satisfied if the number of mesh points that do not satisfy the condition absolute residue < reference absolute residue ε1 is equal to or less than a reference number and the number of mesh points that do not satisfy the condition relative residue < reference relative residue ε2 is equal to or less than a reference number.

[0224] If the convergence condition is satisfied (YES in step S318), the process proceeds to step S321; if the convergence condition is not satisfied (NO in step S318), the process proceeds to step S319.

[0225] The processes in steps S319 and S320 are the same as those in steps S303 and S304 in FIG.

[0226] (Step S321) The gas-liquid interface calculation unit 127 calculates the pressure p, the flow velocity u, the temperature T, the free energy and the solvation free energy μ ｄｉｓ ｋ and the structural information are applied to the Cahn-Hilliard or Allen-Cahn equation to calculate the gas-liquid interface. The structural information is used as a constraint when solving the Cahn-Hilliard or Allen-Cahn equation.

[0227] In the phase-field method, it is assumed that there is a constant thickness between the phases and that physical properties such as density ρ change continuously, and the gas-liquid interface is moved based on the free energy minimization principle, which is a principle that minimizes the free energy possessed by the gas phase, liquid phase, and gas-liquid interface.

[0228] Generally, the phase field method is described by the Allen-Cahn equation as shown in equation (29) when the calculation system is a conservative system, and by the Cahn-Hilliard equation as shown in equation (30) when the calculation system is a non-conservative system. The gas-liquid interface calculation unit 127 may use either equation, but the mobility M φ If there is a high possibility that is location-dependent due to reaction and dissolution, the Cahn-Hilliard equation can be used. The PF variable φ takes values ​​between -1 and 1. The PF variable φ is determined by the volume ratio V of the liquid and gas components as shown in the above equation (22). ｌｉｑ , V ｇａｓ When φ = -1, the liquid phase is 100%. When φ = 1, the gas phase is 100%. When φ = a value other than -1 or 100, the PF variable φ represents the gas-liquid interface.

[0229]

[0230] In equations (29), (30), and (31), u is the flow velocity, F is the free energy of the entire system, and f ｇａｓ＿ｌｉｑ is the free energy per unit volume of the gas-liquid interface, F ｇａｓ＿ｌｉｑ is the free energy of the gas-liquid interface.ｇａｓ ｋ is the free energy per unit volume of component k in the gas phase, f ｌｉｑ ｌ is the free energy per unit volume of component l in the liquid phase. ｇａｓ ｋ is the free energy of component k in the gas phase, F ｌｉｑ ｌ is the free energy of component l in the liquid phase. ｄｉｓ ｋ is the free energy per unit volume that changes when component k in the gas phase dissolves in the liquid phase. M φ is the mobility of the gas-liquid interface. ε ｐｆ is the thickness of the gas-liquid interface.

[0231] Equation (32) is the free energy f per unit volume of the gas-liquid interface in the phase field method. ｇａｓ＿ｌｉｑ can be expressed as the square of the density change. The flow velocity u uses the value of the fluid calculation unit 122.

[0232] Next, the parameter ε included in the phase field method ｐｆ ２ and mobility M φ A method for determining the above will be described.

[0233] First, the parameter ε ｐｆ ２ In the equilibrium state in equation (29) or (30), ∂φ / ∂t = 0, so the left side is 0. In addition, since the chemical potentials of the gas phase and the liquid phase are equal at the gas-liquid interface, equation (33) is established from equations (29), (30), (31), and (32). Therefore, according to equation (33), the parameter ε ｐｆ ２ is determined.

[0234]

[0235] Next, the mobility M φ Assuming that the gas-liquid interface moves at a constant speed in an equilibrium state after a sufficient amount of time has passed, the speed of movement of the gas-liquid interface, u ａｖｇ can be expressed by equation (34) using the mobility parameter α.

[0236]

[0237] The mobility parameter α may be an experimental value, or may generally be calculated with α=1.

[0238] From equations (29) and (30), the following equation holds in the equilibrium state. Therefore, the mobility M φ is derived by equation (35).

[0239]

[0240]

[0241] The gas-liquid interface is determined by simulation using the phase-field method, in which the free energy is described by the sum of many-body potentials based on cluster expansion.

[0242] Conventionally, the free energy f ｌｉｑ , f ｇａｓ Conventionally, an easily calculated function that ignores the physical properties of the target molecule has been used for . Therefore, the present disclosure focuses on the virial expansion, which is a cluster expansion, for the gas-liquid interface. That is, the present disclosure applies the virial expansion to the equation of state of a certain component i using the density ρ and temperature T. Equations (36) and (37) are the equations of state to which the virial expansion has been applied.

[0243]

[0244] The advantage of the virial expansion is that it can describe the equation of state for pressure P using only the potential between N bodies, and in principle, it can describe the equation of state precisely based on physical properties. Another advantage of the virial expansion is that it can describe the intermolecular potential U ｋｌ If it can be expressed as a function, the intermolecular potential U ｋｌ The advantage of this is that the pressure p can be expressed as the sum of m. More preferably, when N=4 or more, the pressure p obtained here is statistically related to the Helmholtz free energy with higher accuracy. Therefore, under constant temperature conditions, m ｇａｓ＿ｉ is the molar ratio of component i in the gas phase, m ｌｉｑ＿ｊ is the molar ratio of component j in the liquid phase, P ｌｉｑ ｉ is the partial pressure of component i in the liquid phase, P ｇａｓｊ If is the partial pressure of component j in the gas phase, then equations (38) and (39) hold.

[0245]

[0246] If an equation of state including pressure p and density ρ is obtained in equations (36) and (37) through experiments, model calculations assuming van der Waals gas, or molecular dynamics calculations, the coefficient β of the virial expansion can be calculated. Ｎ and N of the N-body potential are determined. Then, the coefficients of the virial expansion β Ｎ and N of the N-body potential are determined, the free energies of the liquid phase and the gas phase can be determined in principle.

[0247] Regarding dissolution between liquid and gas phases, we focus on the component k that dissolves from the gas phase into the liquid phase, and define the free energy per unit volume that changes due to dissolution of component k in a certain gas phase as f ｄｉｓ ｋ The solvation free energy of component k is μ ｄｉｓ ｋ When component k dissolves in the liquid phase, the concentration c ｇａｓ＿ｄｉｖ ｋ Since only the free energy of dissolution is involved, we can use this to calculate the free energy of dissolution per unit volume, f ｄｉｓ is written in equation (41).

[0248]

[0249] Solvation free energy μ of component k ｄｉｓ ｋ is determined based on experimental values, model values ​​such as the van der Waals model, or simulations such as molecular dynamics. ｄｉｓ is determined.

[0250] As described above, the gas-liquid interface calculation unit 127 calculates the free energy (f ｇａｓ＿ｌｉｑ , F ｇａｓ＿ｌｉｑ , f ｇａｓ ｋ , f ｌｉｑ ｌ , F ｇａｓ ｋ , F ｌｉｑｌ , f ｄｉｓ ｋ ) is determined, and the parameter ε ｐｆ ２ , the moving speed u of the gas-liquid interface shown in equation (34) ａｖｇ , and the mobility M shown in (35) φ Finally, the PF variable φ indicating the gas-liquid interface is determined from equation (29) or equation (30).

[0251] (Step S322) The gas-liquid interface calculation unit 127 determines whether the calculation result of step S321 satisfies the convergence condition. As in step S302, the gas-liquid interface calculation unit 127 calculates the absolute residual and the relative residual of the nth solution for each mesh point. The solution here is the solution of the Allen-Cahn equation shown in equation (29), i.e., the PF variable φ. Alternatively, the solution here is the solution of the Cahn-Hilliard equation shown in equation (30), i.e., the PF variable φ.

[0252] The gas-liquid interface calculation unit 127 may determine that the convergence condition is satisfied if the absolute residue is less than the reference absolute residue ε1 and the relative residue is less than the reference relative residue ε2. Values ​​set in advance by the user are adopted as the reference absolute residue ε1 and the reference relative residue ε2. Values ​​suitable for the solution of the gas-liquid interface calculation unit 127 are adopted as the reference absolute residue ε1 and the reference relative residue ε2. The gas-liquid interface calculation unit 127 may determine that the convergence condition is satisfied if the number of mesh points that do not satisfy the condition absolute residue < reference absolute residue ε1 is less than or equal to a reference number and the number of mesh points that do not satisfy the condition relative residue < reference relative residue ε2 is less than or equal to a reference number.

[0253] If the convergence condition is satisfied (YES in step S322), the process proceeds to step S325; if the convergence condition is not satisfied (NO in step S322), the process proceeds to step S323.

[0254] The processes in steps S323 and S324 are the same as those in steps S303 and S304 in FIG.

[0255] (Step S325) The overall calculation unit 128 calculates, for each mesh point, the step absolute residual and the step relative residual of the solutions in the fluid calculation unit 122, the chemical reaction calculation unit 123, the heat transfer calculation unit 124, the diffusion calculation unit 125, the physical property calculation unit 126, and the gas-liquid interface calculation unit 127. The definitions of the step absolute residual and the step relative residual are the same as those of the absolute residual and the relative residual described in step S302. After step S325 is completed, the process proceeds to step S4 in FIG. 3, where it is determined whether the step absolute residual and the step relative residual satisfy the step convergence condition.

[0256] (Method of solving equations) Either an explicit method or an implicit method may be adopted as a numerical analysis method for solving equations in each of the fluid calculation unit 122, the chemical reaction calculation unit 123, the heat transfer calculation unit 124, the diffusion calculation unit 125, the physical property calculation unit 126, and the gas-liquid interface calculation unit 127. An explicit method is adopted when chemical changes, dissolution, and changes in the gas-liquid interface change dramatically in a relatively short period of time. Examples of explicit methods include the Runge-Kutta method and the Adams-Bashforth method.

[0257] On the other hand, when chemical changes, dissolution, and changes in the gas-liquid interface do not change significantly in a short time but do change over a long time, an implicit method is used. Examples of implicit methods include BDF (registered trademark) (Backward Differentiation Formula) and the generalized α method. BDF or the generalized α method, which have high numerical stability as a solution method, are more desirable.

[0258] (Convergence Parameter) In the present disclosure, the Newton method or the quasi-Newton method is adopted as a solution convergence method. As the quasi-Newton method, generally, the DFP (Davidon-Fletcher-Powell) method, the BFGS (Broyden-Fletcher-Goldfarb-Shanno) method, the L-BFGS method, the SR1 (Symmetric Rank One) method, etc. can be adopted. In the Newton method, when a solution of f(x) = 0 is obtained, the nth solution is x = x ｎ Then, it is expressed by equation (42).

[0259]

[0260] In the Newton method, for example, a damping factor λ is used as a parameter for the convergence of the solution. The damping factor λ is shown in equation (43).

[0261]

[0262] Equation (43) is an equation incorporating a damping factor λ into equation (42). By using the damping factor λ, divergence of the solution is prevented, and the solution can be obtained in a shorter time. The convergence parameter is, for example, the Jacobian matrix f'(x ｎ ) may be the number of updates of the Jacobian matrix f'(x ｎ ) increases the calculation cost, but improves the convergence. Anderson acceleration may also be used to find the convergence parameter.

[0263] As described above, according to the information processing system 1 of the present embodiment, the gas-liquid interface is calculated taking into consideration the first reaction property value calculated based on the reaction information, so that the gas-liquid interface can be predicted accurately. Furthermore, since the reaction includes dissolution and the gas-liquid interface is calculated so that the difference between the free energy of the gas phase and the free energy of the liquid phase is small, the gas-liquid interface can be predicted more accurately.

[0264] (Example) In this example, the present disclosure is applied to the combustion of methane in an engine. In this example, polyolefin is used as the engine oil. The combustion of methane was calculated using the following formula:

[0265]

[0266] In this example, methane was used as the gas phase and a polyolefin with a molecular weight of 2000 was used as the liquid phase. The initial volume ratio of methane, oxygen, and polyolefin charged into the engine was set to 1:0.5:0.5. The initial temperature of the methane, oxygen, and polyolefin was set to 300 K. A cylindrical flow path with a radius of 1 cm and a height of 10 cm was used. Only the polyolefin was flowed from the inlet at a temperature of 300 K and a velocity of 5 cm / s. A spherical surface (heat source) with a radius of 0.2 cm was installed midway through the flow path. The temperature of this spherical surface was set to 400 K after 0.2 seconds. Using this temperature setting, gaseous methane was combusted. The calculations assumed that all water and carbon dioxide produced by combustion were gaseous. Because a portion of the methane dissolves in the polyolefin, the chemical potential was calculated in advance using molecular dynamics calculations. Equation (41) was determined using this chemical potential.

[0267] Figure 16 shows the state of the gas-liquid interface after 0.5 seconds, based on calculation results using equations (3), (4), (5), (6), (7), (8), (9), (10), (11), (12), (13), (14), (15), (16), (17), (27), (33), (34), (35), (36), (37), (38), (39), (40), and (41). Equation (1) incorporates the effect of flow due to surface tension at the gas-liquid interface. Equation (1) prevents tiny bubbles from being crushed and disappearing. Equation (6) shows the change in chemical potential at the gas-liquid interface. Equation (6) determines the complex shape of the gas-liquid interface caused by reactions, etc. Equation (12) incorporates the effect of reaction heat generated by chemical reactions. Equation (13) incorporates the effect of methane dissolving in polyolefin, increasing the liquid phase. Equation (14) incorporates the effects of reaction heat and dissolution heat on the gas-liquid interface, and the temperature is determined from the heat conduction equation. Equations (27), (33), and (34) incorporate the free energies of methane, oxygen, and polyolefin, and a precise gas-liquid interface can be obtained by expressing the gas-liquid interface in principle. Furthermore, by incorporating equation (36), the dissolution effect of methane in polyolefin can be incorporated.

[0268] (Comparative Example) In the comparative example, the same calculations as in the example were performed without including formulas (3), (4), (5), (6), (7), (8), (9), (10), (11), (12), (13), (14), (15), (16), (17), (27), (33), (34), (35), (36), (37), (38), (39), (40), and (41). In formula (14), the calculation does not include the heat due to the combustion reaction. Therefore, in the comparative example, the calorific value R ｔｏｔ = 0, which is a heat conduction equation that does not depend on the phase field method. In equation (22), since the dissolution effect is not included at all, the calculation was performed with the right side = 0. Figure 17 shows the simulation results of the comparative example.

[0269] Comparing the simulation results after 0.5 seconds between the Example and the Comparative Example, the Example had a large growth of gas phase 1007 as shown in Fig. 16. In contrast, in the Comparative Example, the gas phase 1007 had almost completely disappeared, and the flow path was filled with liquid phase 1008 as shown in Fig. 17. From these calculation results, it was confirmed that the Example accurately reproduced a gas-liquid two-layer flow incorporating the reaction.

[0270] The present disclosure is useful in the technical field of simulating gas-liquid multiphase flow.

Claims

1. An information processing method performed by a computer to identify a gas-liquid interface of a gas-liquid multiphase flow flowing in a flow path, the information processing method comprising: acquiring first gas-liquid multiphase flow information including the free energy of each of one or more substances included in the gas-liquid multiphase flow; acquiring reaction information related to reactions occurring between one or more substances included in the gas phase and one or more substances included in the liquid phase; identifying the gas-liquid interface of the gas-liquid multiphase flow using the free energy of the gas phase and the free energy of the liquid phase; and outputting gas-liquid interface information indicating the identified gas-liquid interface.

2. The information processing method according to claim 1, further comprising: acquiring structural information including a structure including the flow path and boundary conditions for calculating the gas-liquid interface of the gas-liquid multiphase flow; and calculating one or more reaction property values ​​obtained by the reaction based on the first gas-liquid multiphase flow information, the structural information, and the reaction information.

3. The information processing method of claim 2, wherein identifying the gas-liquid interface includes calculating the gas-liquid interface by executing a simulation using a phase field method with the structural information, the reaction information, and the one or more reaction property values ​​as inputs, and wherein the simulation using the phase field method describes the free energy as a sum of many-body potentials based on a cluster expansion.

4. An information processing method according to claim 1 or 2, wherein the reaction includes dissolving the one or more substances contained in the gas phase into the one or more substances in the liquid phase, and identifying the gas-liquid interface includes identifying the gas-liquid interface so that the difference between the free energy of the gas phase and the free energy of the liquid phase is small.

5. The information processing method according to claim 2, wherein identifying the gas-liquid interface includes calculating surface tension based on the temperature or dissolution calculated based on the reaction information, and calculating the gas-liquid interface so that the surface tension matches the external pressure.

6. The information processing method of claim 2, wherein calculating the one or more reaction property values ​​further comprises: calculating at least one of the amount of reaction component and the heat of reaction obtained by the reaction based on the first gas-liquid multiphase flow information, the structural information, and the reaction information; and calculating the temperature of the gas-liquid multiphase flow due to the reaction based on the first gas-liquid multiphase flow information, the structural information, and at least one of the amount of reaction component and the heat of reaction.

7. An information processing method as described in claim 2, further comprising calculating second gas-liquid multiphase flow information including one or more second property values ​​of the gas-liquid multiphase flow based on the first gas-liquid multiphase flow information, the calculated one or more reactant property values, and the identified gas-liquid interface, wherein calculating the one or more reactant property values ​​comprises calculating the one or more reactant property values ​​obtained by the reaction based on the second gas-liquid multiphase flow information, the structural information, and the reaction information, and calculating the gas-liquid interface of the gas-liquid multiphase flow comprises calculating the gas-liquid interface of the gas-liquid multiphase flow based on the second gas-liquid multiphase flow information, the structural information, the reaction information, and the one or more reactant property values.

8. The information processing method according to claim 1 or 2, wherein the output gas-liquid interface information indicates a gas-liquid interface calculated at a predetermined elapsed time.

9. An information processing method according to claim 2, wherein acquiring the structural information includes acquiring a plurality of mesh points included in the structure, and further includes outputting second gas-liquid multiphase flow information including one or more second physical property values ​​of the gas-liquid multiphase flow calculated for each of the plurality of mesh points.

10. An information processing method according to claim 1 or 2, further comprising outputting second gas-liquid multiphase flow information including one or more second physical property values ​​of the gas-liquid multiphase flow calculated at a predetermined elapsed time in the simulation.

11. The information processing method according to claim 2, wherein identifying the gas-liquid interface includes calculating the first gas-liquid multiphase flow information, the structural information, the reaction information, and the one or more reaction property values ​​by substituting them into a Cahn-Hilliard or Allen-Cahn equation.

12. An information processing method according to claim 2, wherein acquiring the structural information includes acquiring a plurality of mesh points included in the structure, and further includes determining whether or not a calculated value of the gas-liquid interface of the gas-liquid multiphase flow calculated for each of the plurality of mesh points is below a predetermined reference value, and changing a convergence parameter used in calculating the gas-liquid interface of the gas-liquid multiphase flow in accordance with the determination result, and identifying the gas-liquid interface includes identifying the gas-liquid interface of the gas-liquid multiphase flow based on the changed convergence parameter.

13. A system for providing physical property information of a gas-liquid multiphase flow flowing through a flow path, comprising: a display control unit that causes a display unit to display a first image including a flow path model representing the flow path and data representing mesh points present in the flow path model or time information indicating a predetermined elapsed time; and an image generation unit that, upon receiving a selection of the data representing the mesh points or the time information, generates a second image including data representing calculated values ​​of one or more physical properties of the gas-liquid multiphase flow.

14. A physical property information provision system as described in claim 13, wherein, when the image generation unit receives a selection of data indicating the mesh point, it generates the second image including data representing calculated values ​​of one or more physical properties of the gas-liquid multiphase flow at the selected mesh point.

15. An information processing program for identifying a gas-liquid interface of a gas-liquid multiphase flow flowing in a flow path, the information processing program causing a computer to execute the following operations: acquire first gas-liquid multiphase flow information including the free energy of each of one or more substances included in the gas-liquid multiphase flow; acquire reaction information related to reactions occurring between one or more substances included in the gas phase and one or more substances included in the liquid phase; identify the gas-liquid interface of the gas-liquid multiphase flow using the free energy of the gas phase and the free energy of the liquid phase; and output gas-liquid interface information indicating the identified gas-liquid interface.

Images