WO2021246378A1 - シミュレーション装置、シミュレーション方法、プログラム - Google Patents

シミュレーション装置、シミュレーション方法、プログラム Download PDF

Info

Publication number
WO2021246378A1
WO2021246378A1 PCT/JP2021/020745 JP2021020745W WO2021246378A1 WO 2021246378 A1 WO2021246378 A1 WO 2021246378A1 JP 2021020745 W JP2021020745 W JP 2021020745W WO 2021246378 A1 WO2021246378 A1 WO 2021246378A1
Authority
WO
WIPO (PCT)
Prior art keywords
coarse
parameter
particles
grained
powder
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Ceased
Application number
PCT/JP2021/020745
Other languages
English (en)
French (fr)
Inventor
元彬 猿渡
英也 仲村
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Sumitomo Metal Mining Co Ltd
University Public Corporation Osaka
Original Assignee
Sumitomo Metal Mining Co Ltd
University Public Corporation Osaka
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Priority claimed from JP2020162719A external-priority patent/JP7101387B2/ja
Application filed by Sumitomo Metal Mining Co Ltd, University Public Corporation Osaka filed Critical Sumitomo Metal Mining Co Ltd
Priority to EP21818068.5A priority Critical patent/EP4160184B1/en
Priority to CN202180039369.1A priority patent/CN115698672B/zh
Priority to KR1020227041975A priority patent/KR102660215B1/ko
Publication of WO2021246378A1 publication Critical patent/WO2021246378A1/ja
Priority to US18/058,808 priority patent/US11892388B2/en
Anticipated expiration legal-status Critical
Ceased legal-status Critical Current

Links

Images

Classifications

    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16C—COMPUTATIONAL CHEMISTRY; CHEMOINFORMATICS; COMPUTATIONAL MATERIALS SCIENCE
    • G16C60/00—Computational materials science, i.e. ICT specially adapted for investigating the physical or chemical properties of materials or phenomena associated with their design, synthesis, processing, characterisation or utilisation
    • G—PHYSICS
    • G06—COMPUTING OR CALCULATING; COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/25—Design optimisation, verification or simulation using particle-based methods
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N15/00—Investigating characteristics of particles; Investigating permeability, pore-volume or surface-area of porous materials
    • G01N15/10—Investigating individual particles
    • G01N15/1012—Calibrating particle analysers; References therefor
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N15/00—Investigating characteristics of particles; Investigating permeability, pore-volume or surface-area of porous materials
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N15/00—Investigating characteristics of particles; Investigating permeability, pore-volume or surface-area of porous materials
    • G01N15/10—Investigating individual particles
    • G01N15/14—Optical investigation techniques, e.g. flow cytometry
    • G01N15/1425—Optical investigation techniques, e.g. flow cytometry using an analyser being characterised by its control arrangement
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16C—COMPUTATIONAL CHEMISTRY; CHEMOINFORMATICS; COMPUTATIONAL MATERIALS SCIENCE
    • G16C20/00—Chemoinformatics, i.e. ICT specially adapted for the handling of physicochemical or structural data of chemical particles, elements, compounds or mixtures
    • G16C20/30—Prediction of properties of chemical compounds, compositions or mixtures
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16C—COMPUTATIONAL CHEMISTRY; CHEMOINFORMATICS; COMPUTATIONAL MATERIALS SCIENCE
    • G16C20/00—Chemoinformatics, i.e. ICT specially adapted for the handling of physicochemical or structural data of chemical particles, elements, compounds or mixtures
    • G16C20/60—In silico combinatorial chemistry
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16C—COMPUTATIONAL CHEMISTRY; CHEMOINFORMATICS; COMPUTATIONAL MATERIALS SCIENCE
    • G16C20/00—Chemoinformatics, i.e. ICT specially adapted for the handling of physicochemical or structural data of chemical particles, elements, compounds or mixtures
    • G16C20/90—Programming languages; Computing architectures; Database systems; Data warehousing
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N15/00—Investigating characteristics of particles; Investigating permeability, pore-volume or surface-area of porous materials
    • G01N2015/0003—Determining electric mobility, velocity profile, average speed or velocity of a plurality of particles
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N15/00—Investigating characteristics of particles; Investigating permeability, pore-volume or surface-area of porous materials
    • G01N15/10—Investigating individual particles
    • G01N15/1012—Calibrating particle analysers; References therefor
    • G01N2015/1016—Particle flow simulating, e.g. liquid crystal cell
    • G—PHYSICS
    • G06—COMPUTING OR CALCULATING; COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/10—Numerical modelling
    • G—PHYSICS
    • G06—COMPUTING OR CALCULATING; COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/08—Thermal analysis or thermal optimisation
    • G—PHYSICS
    • G06—COMPUTING OR CALCULATING; COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces

Definitions

  • the present invention relates to a simulation apparatus, a simulation method, and a program.
  • Patent Document 1 describes an input device for inputting simulation conditions and an input device.
  • An output device that outputs simulation results and It has a processing device for analyzing the behavior of powder or granular material containing a plurality of particles having different sizes based on the simulation conditions input from the input device.
  • the processing device is The coarse-grained powder is based on the value of the parameter that defines the particle size distribution of the powder or granular material to be simulated and the value of the coarse-graining coefficient that is the reference for coarse-graining the particles, which is input from the input device.
  • a simulation device for powder or granular materials is disclosed as disclosed in a simulation device or the like that associates the behavior of particles obtained by simulation with the input coarse-graining coefficient value and outputs the result to the output device. ing.
  • Discrete Element Method is used to analyze the behavior of powders (powder and granular materials) containing multiple particles for the purpose of improving processes in factories and reducing the number of test man-hours when considering manufacturing processes. ) It has been conventionally done by calculation or the like.
  • Discrete element method calculation is a simulation technique that describes the motion of the entire powder by solving the equation of motion for each particle.
  • a simulation device using a coarse-graining method in which a particle group composed of a plurality of particles is used as one coarse-grained particle has been conventionally studied (see Patent Document 1).
  • a simulation device using the coarse-graining method it is required to appropriately select the parameters related to the coarse-grained particles used in the calculation in order to obtain accurate analysis results. Therefore, there has been a demand for a new simulation device capable of selecting and setting parameters related to coarse-grained particles by a new method and analyzing the behavior of powder containing a plurality of particles.
  • one aspect of the present invention is to provide a new simulation apparatus capable of analyzing the behavior of a powder containing a plurality of particles.
  • a simulation device for analyzing the behavior of powder containing multiple particles A first parameter acquisition unit that acquires a first parameter including parameters related to the powder, and a first parameter acquisition unit.
  • a second parameter calculation unit that calculates a second parameter that is a parameter for the coarse-grained particles when a particle group composed of a plurality of the particles is coarse-grained and made into one coarse-grained particle.
  • It has a coarse-grained particle behavior analysis unit that analyzes the behavior of the coarse-grained particles based on the first parameter and the second parameter.
  • the second parameter calculation unit provides a simulation device for calculating the second parameter by using a solution of a characteristic equation using the relationship between the elastic energy of the particle group and the elastic energy of the coarse-grained particles. ..
  • the present embodiment Specific examples of the simulation apparatus, simulation method, and program according to the embodiment of the present disclosure (hereinafter referred to as “the present embodiment”) will be described below with reference to the drawings. It should be noted that the present invention is not limited to these examples, and is indicated by the scope of claims, and is intended to include all modifications within the meaning and scope equivalent to the scope of claims. 1.
  • Simulation device (1) Particle coarse-graining and parameters used to calculate the particle behavior of coarse-grained particles (1-1) Particle coarse-graining
  • this embodiment A method of coarse-graining a particle group composed of a plurality of particles and a method of calculating parameters related to coarse-grained particles, which are coarse-grained particles, that can be used in the morphological simulation device will be described below.
  • the calculation load increases as the number of particles handled increases. For this reason, when analyzing the behavior of powder on a large scale such as a plant used in a factory, the amount of calculation becomes enormous, and it becomes difficult to perform the calculation in reality.
  • one particle group 11 composed of a plurality of particles shown in FIG. 1A is as shown in FIG. 2A.
  • the particle group 11 composed of a plurality of particles, the cubic, longitudinal, lateral, two by two in the height direction, and are arranged 2 3 in total.
  • the number of particles arranged in the direction of one side, that is, 2 is defined as the coarse-grained magnification.
  • FIG. 1A is a side view of the state when the particle swarm 11 collides with the wall surface 12.
  • the magnitude of the force applied to the particle group 11 can be expressed by the following equation (1).
  • ⁇ in the equation (1) is the coarse-grained magnification, and means the number of particles arranged in one side when the particle group 11 is regarded as one coarse-grained particle as described above.
  • ⁇ 2.
  • m means the mass of each particle 11A and 11B
  • a G means the acceleration of the center of gravity of the particle group 11
  • ⁇ w means the coefficient of restitution calculated from the restitution coefficient between the wall surface 12 or the external particle and the particle 11A.
  • Equation (A) The relationship between the coefficient of restitution e and the coefficient of viscosity ⁇ used when calculating the above-mentioned ⁇ w or the like from the coefficient of restitution can be expressed by the following equation (A).
  • m * means reduced mass and K means spring constant.
  • the particle group 11 composed of the eight particles shown in FIG. 1A is designated as one coarse-grained particle 21.
  • the force received by the coarse-grained particles 21 can be expressed by the following equation (2).
  • F cw in the formula (2) is the force that the coarse-grained particles 21 receive from the wall surface 12 or the external particles as shown in FIG. 2B, and ⁇ cw is the overlap with the wall surface 12 or the external particles of the coarse-grained particles 21.
  • the amount and ⁇ cw mean the viscosity coefficient calculated from the repulsion coefficient between the coarse-grained particles 21 and the wall surface 12 or the external particles.
  • FIG. 2B is a side view of the state of the coarse-grained particles 21 when they collide with the wall surface 12.
  • the corresponding parameter is obtained from the above-mentioned equation (1) in which the calculation is performed for the particle group 11 and the above-mentioned equation (2) in which the calculation is performed for the coarse-grained particle 21 in which the particle group is coarse-grained.
  • equations (3) and (4) indicating that they match are derived.
  • the force applied to each particle can be expressed as the following equations (5) to (7) by using the particle overlap amounts ⁇ w , ⁇ p , and ⁇ cw. can.
  • K w is the spring coefficient between the particle 11A and the wall surface 12 or the external particle
  • K p is the spring coefficient of the internal particle of the particle group 11
  • K cw is the coarse-grained particle 21. It means the spring coefficient with the wall surface 12 or the external particles, respectively.
  • the overlap amount [delta] CW of the wall 12 or the external particles of the coarse-grained particles 21 the overlap amount of the particles constituting the coarse-grained front of the particle group 11 calculated I also know that I can do it.
  • K r can be calculated by the characteristic equation using the elastic energy of the coarse-grained front of the particle group 11 in a collision, the relationship between the elastic energy of the coarse-grained particles 21. Specifically, for example, assuming the total elastic energy coarse-grained front of the particle group 11, the elastic energy of the entire coarse-grained particles are equal, it can be calculated K r to create a characteristic equation.
  • the elastic energy of the particle group 11 and the elastic energy of the coarse-grained particles during collision with the wall surface 12 exert a force applied to the particles 11A and 11B constituting the above-mentioned particle group 11 and a force applied to the coarse-grained particles. It can be calculated by integrating the expressed equations (5) to (7) with the distance of the overlap amount.
  • the following equation (12) can be obtained by using the elastic energy of the entire particle group 11 before coarse graining and the elastic energy of the entire coarse-grained particles.
  • the above equation (12) can be transformed into the following equation (13) by using the above-mentioned equations (8) to (11).
  • Equation (13) is a characteristic equation of Kr in the vertical direction.
  • K r is a parameter related to the overlap amount of the particles constituting the coarse-grained particles, and a storm of previous particle group 11, the coarse-grained particles It is a parameter that governs the behavior. Therefore, by keeping seeking K r in advance by a characteristic equation, such as to calculate the overlap amount of the coarse-grained before particles from the overlap amount of coarse-grained particles after the parameters for coarse-grained particles It is possible to calculate and calculate the behavior of coarse-grained particles.
  • equation (15) when Kr is set, ⁇ w and ⁇ p can be expressed as equations (16) and (17), and the particle group before coarse graining.
  • equation (18) is obtained.
  • equation (19) which is a characteristic equation in the tangential direction, can be obtained.
  • a linear spring model was used for the tangential contact model.
  • the elastic energy calculation formula differs depending on the contact model, but the elastic energy can be appropriately calculated by changing the characteristic equation as necessary.
  • the heat transfer of particles is expressed by the following equations (20) and (21) using thermal conductivity.
  • Equation (20) Q is the heat flow
  • h is the heat flow coefficient
  • [Delta] T is the temperature difference between the wall surface or particle
  • k p is the thermal conductivity of the internal particle group 11, a means contact radius of the particle group 11, respectively.
  • the contact radius depends on the amount of overlap of each particle. Then, as described above, the amount of overlap before and after coarse graining is related to the solution of the above-mentioned vertical characteristic equation. Therefore, vertical characteristic equation solutions K following equations r with rewriting the heat transfer equation of (22), so that equation (23).
  • the heat transfer equation of the coarse-grained particles can be expressed by the following equations (24) and (25).
  • Q c is the heat flow rate of the coarse-grained particles
  • h c is the heat flow coefficient of the coarse-grained particles
  • ⁇ T c is the coarse-grained particles 21 and the wall surface 12 or external particles.
  • the temperature difference of k w ′ means the thermal conductivity between the wall surface 12 and the coarse-grained particles 21, and k p ′ means the thermal conductivity between the external particles and the coarse-grained particles.
  • the simulation device of the present embodiment is a simulation device for analyzing the behavior of a powder containing a plurality of particles, and can have the following members.
  • a first parameter acquisition unit that acquires a first parameter including parameters related to powder.
  • a second parameter calculation unit that calculates a second parameter that is a parameter for coarse-grained particles when a particle group composed of a plurality of particles is coarse-grained and made into one coarse-grained particle.
  • a coarse-grained particle behavior analysis unit that analyzes the behavior of coarse-grained particles based on the first parameter and the second parameter.
  • the second parameter calculation unit calculates the second parameter using the solution of the characteristic equation using the relationship between the elastic energy of the particle group and the elastic energy of the coarse-grained particles.
  • the simulation device 30 of the present embodiment is composed of, for example, an information processing device (computer), and is physically a CPU (Central Processing Unit) which is an arithmetic processing unit. : Processor) 31, including RAM (Random Access Memory) 32 and ROM (Read Only Memory) 33 which are main storage devices, auxiliary storage device 34, input / output interface 35, display device 36 which is an output device, and the like. It can be configured as a computer system. These are interconnected by a bus 37.
  • the auxiliary storage device 34 and the display device 36 may be provided externally.
  • the CPU 31 controls the overall operation of the simulation device 30 and performs various information processing.
  • the CPU 31 executes a simulation method or a program (simulation program) stored in the ROM 33 or the auxiliary storage device 34, which will be described later, to calculate a second parameter which is a parameter for coarse-grained particles, and to perform coarse-graining. It is possible to analyze the behavior of particles.
  • the RAM 32 may include a non-volatile RAM that is used as the work area of the CPU 31 and stores major control parameters and information.
  • the ROM 33 can store a program (simulation program) or the like.
  • the auxiliary storage device 34 is a storage device such as an SSD (Solid State Drive) or an HDD (Hard Disk Drive), and can store various data, files, and the like necessary for the operation of the simulation device.
  • SSD Solid State Drive
  • HDD Hard Disk Drive
  • the input / output interface 35 includes both a user interface such as a touch panel, a keyboard, a display screen, and operation buttons, and a communication interface that takes in information from an external data recording server and outputs analysis information to other electronic devices. ..
  • the display device 36 is a monitor display or the like.
  • the display device 36 displays an analysis screen, and the screen is updated according to the input / output operation via the input / output interface 35.
  • Each function of the simulation device 30 shown in FIG. 3 is performed by reading a program (simulation program) or the like from a main storage device such as a RAM 32 or a ROM 33 or an auxiliary storage device 34 and executing the program by the CPU 31 to obtain data in the RAM 32 or the like. This can be achieved by reading and writing and operating the input / output interface 35 and the display device 36.
  • a program simulation program
  • main storage device such as a RAM 32 or a ROM 33 or an auxiliary storage device 34
  • FIG. 4 shows a functional block diagram of the simulation device 30 of the present embodiment.
  • the simulation device 30 can have a reception unit 41, a processing device 42, and an output unit 43.
  • Each of these parts is an information processing device such as a personal computer provided with a CPU, a storage device, various interfaces, etc. of the simulation device 30, and the CPU is stored in advance, for example, by executing a simulation method or a program described later.
  • Software and hardware work together.
  • the configuration of each part will be described below.
  • A Reception unit
  • the reception unit 41 receives input of commands and data from a user related to processing executed by the processing device 42.
  • the reception unit 41 includes a keyboard and mouse for users to operate and input commands, a communication device for input via a network, a reading device for input from various storage media such as a CD-ROM and a DVD-ROM, and the like. Can be mentioned.
  • B Processing device
  • the processing device 42 can have a first parameter acquisition unit 421, a second parameter calculation unit 422, and a coarse-grained particle behavior analysis unit 423.
  • the processing device may further have an arbitrary member, if necessary, and may also have, for example, an initial setting unit and the like.
  • the first parameter acquisition unit 421 can acquire, for example, a first parameter including a parameter related to the powder to be analyzed.
  • the first parameter may include various parameters required for analysis in addition to the parameters related to the powder. Since the first parameter can be selected according to the content of the analysis (simulation), the specific type thereof is not particularly limited.
  • the first parameter includes various parameters required for the discrete element method calculation. Specifically, for example, the particle size, the number of particles, the Young's ratio, the calculated Time step, the Poisson ratio, the coefficient of friction with the wall surface, and the distance between particles. One or more types selected from the coefficient of friction, the coefficient of rolling friction, the density, and the like can be mentioned.
  • the first parameter may be data recorded in a database or the like, or may be an experimental value obtained by conducting an experiment in advance. Further, the first parameter may be a calculated value calculated by fitting from the experimental result by simulation or the like.
  • B-2 Second Parameter Calculation Unit
  • Coarse-grained particles can be coarse-grained and the number of particles can be reduced. However, the coarse-grained particles differ from the individual particles constituting the particle group before coarse-graining in various parameters such as mass. Therefore, it is necessary to calculate and set the parameters required for the calculation of the coarse-grained particles.
  • the K r is a solution of the derived characteristic equation by using the relationship between the elastic energy of the coarse-grained particles, it calculates a second parameter.
  • the above-mentioned equation Kr characteristic equation in the vertical direction
  • K r is a parameter that governs the behavior of coarse-grained particles
  • the use of K r can be calculated various parameters related to the behavior of the coarse-grained particles.
  • the type of the second parameter used in the coarse-grained particle behavior analysis unit described later is not particularly limited because it can be selected according to the content of the analysis.
  • the second parameter can also include the thermal conductivity of the coarse-grained particles.
  • second parameter calculating unit uses the solution K r of the characteristic equation described above can calculate the thermal conductivity.
  • (B-3) Coarse-grained particle behavior analysis unit
  • the coarse-grained particle behavior analysis unit 423 uses the first parameter acquired by the first parameter acquisition unit 421 and the second parameter calculated by the second parameter calculation unit 422. Therefore, the behavior of coarse-grained particles can be analyzed. Specifically, it can be calculated using the discrete element method and the behavior of coarse-grained particles can be analyzed. By analyzing the behavior of coarse-grained particles, the behavior of powder can be analyzed.
  • the behavior referred to here includes not only a change in position due to the movement of coarse-grained particles but also a change in state such as a temperature change.
  • B-4) Initial setting unit The initial setting unit (not shown) initializes the positions of the particles that make up the powder to be analyzed, and also determines the analysis conditions, for example, the temperature of the region where the powder is placed, if necessary. Etc. can be set. For example, when the initial conditions are set in advance in the program or the like used when the coarse-grained particle behavior analysis unit 423 analyzes the behavior of the coarse-grained particles, or when the first parameter acquisition unit 421 acquires the particles. , The initial setting unit may not be provided.
  • (C) Output unit The output unit 43 may have a display or the like.
  • the simulation result obtained by the coarse-grained particle behavior analysis unit 423 can be output to the output unit 43.
  • the content of the simulation result to be output is not particularly limited, but for example, the positions of the coarse-grained particles can be output and displayed as an image in time series on the output unit 43. Further, for example, the time-series change of the temperature distribution of the powder can be output and displayed as an image on the output unit 43.
  • the behavior of a powder containing a plurality of particles can be simulated, and its use and the like are not particularly limited.
  • it can be suitably used for simulating the behavior of powder in a rotating body such as a kiln. That is, according to the simulation apparatus of this embodiment, it is possible to analyze the behavior of the powder in the rotating body.
  • the amount of calculation can be suppressed by forming a particle group composed of a plurality of particles into one coarse-grained particle. Therefore, it is possible to suppress the amount of calculation and efficiently perform the calculation even for the behavior of the powder on a large scale such as a plant used in a factory.
  • the simulation device of the present embodiment may further include a powder supply device, a reaction furnace, a control device, and the like in order to carry out various manufacturing processes using powder using the simulation results.
  • the powder supply device examples include a device that can store and discharge powder such as a hopper. Since the powder supply device controls the amount of powder discharged and supplied from the hopper or the like to the reaction furnace, it may further have a supply amount adjusting device such as a feeder or a valve. Based on the simulation results, it is preferable to have a plurality of powder supply devices each containing powders having different physical properties so that powders having desired physical properties, for example, a desired average particle size can be supplied.
  • reaction furnace examples include various reaction furnaces such as a heating furnace, and a rotary furnace, for example, a kiln and the like.
  • the powder supply device and the reactor can be connected by piping.
  • the control device supplies the powder supply device with the desired physical properties, for example, the powder having the desired average particle size. Can be controlled. Further, the control device can control the heating conditions of the reactor based on the result of the behavior of the powder in the obtained reactor. Examples of the heating conditions include temperature conditions in the reactor, atmospheric conditions, heating time, and the like.
  • the powder supply device and the reactor are provided with various sensors to detect the amount of powder supplied, the temperature, etc., so that the measurement results measured at arbitrary timings can be supplied to the control device. You may. In this case, the control device can also control each device based on the obtained measurement data.
  • the control device can have a CPU, a main storage device, an auxiliary storage device, an input / output interface, etc. so that data processing such as control conditions can be performed and communication with a powder supply device and a reactor can be performed.
  • the main storage device include RAM and ROM
  • examples of the auxiliary storage device include SSD and HDD.
  • the input / output interface include a powder supply device and a communication interface for exchanging control signals and data with the reactor.
  • the type of communication interface is not particularly limited. Both wired and wireless communication methods can be used, and examples thereof include a wired LAN (Local Area Network) and a wireless LAN.
  • the powder having the desired physical properties is supplied from the powder supply device, and further, under predetermined heating conditions in the reaction furnace.
  • the reaction ratio of the powder can be increased.
  • the heating conditions and the like can be optimized, the amount of energy used in the reaction can be optimized and the productivity can be improved.
  • the simulation device of the present embodiment has the powder supply device or the like
  • the simulation device can also be referred to as a reaction device or the like.
  • the powder supply device, reactor, and control device can be configured so that they can be separated from the device for analyzing the behavior of powder containing multiple particles, and after the simulation results are reflected in the control device. , May be used separately.
  • Simulation method Next, the simulation method of this embodiment will be described. The simulation method of this embodiment can be carried out using, for example, the simulation device described above. Therefore, some of the matters already described will be omitted.
  • the simulation method of this embodiment relates to a simulation method for analyzing the behavior of a powder containing a plurality of particles.
  • the simulation method of this embodiment can be carried out according to the flow chart shown in FIG. 5, and can have the following steps.
  • the second parameter can be calculated by using the solution of the characteristic equation using the relationship between the elastic energy of the particle group and the elastic energy of the coarse-grained particles.
  • First parameter acquisition step (S1) In the first parameter acquisition step (S1), the first parameter including the parameter related to the powder to be analyzed can be acquired.
  • the first parameter acquisition unit 421 can carry out the first parameter acquisition step.
  • the first parameter can be selected according to the content of the analysis, its specific type is not particularly limited.
  • the first parameter various parameters required for the discrete element method calculation can be mentioned. Since a specific example of the first parameter has already been described in the simulation apparatus, the description thereof will be omitted here.
  • the first parameter may be data recorded in a database or the like, or may be an experimental value obtained by conducting an experiment in advance. Further, the first parameter may be a calculated value calculated by fitting from the experimental result by simulation or the like.
  • Second parameter calculation step (S2) In the simulation method of the present embodiment, in order to suppress the calculation amount, the particle group composed of a plurality of particles contained in the powder is coarse-grained into one coarse-grained particle, and the number of particles is reduced. You can make calculations.
  • the K r is a solution of the derived characteristic equation by using the relationship between the elastic energy of the coarse-grained particles, it calculates a second parameter.
  • the above-mentioned equation Kr characteristic equation in the vertical direction
  • K r is a parameter that governs the behavior of coarse-grained particles
  • the use of K r can be calculated various parameters related to the behavior of the coarse-grained particles.
  • the second parameter calculation unit 422 can carry out the second parameter calculation step.
  • the type of the second parameter used in the coarse-grained particle behavior analysis step described later is not particularly limited because it can be selected according to the content of the analysis.
  • the second parameter can also include the thermal conductivity of the coarse-grained particles.
  • second parameter calculating step using a solution K r of the characteristic equation described above can calculate the thermal conductivity.
  • (3) Coarse-grained particle behavior analysis step (S3) In the coarse-grained particle behavior analysis step (S3), the coarse-grained particles are used by using the first parameter obtained in the first parameter acquisition step (S1) and the second parameter calculated in the second parameter calculation step (S2). Behavior can be analyzed. Specifically, it can be calculated using the discrete element method and the behavior of coarse-grained particles can be analyzed. By analyzing the behavior of coarse-grained particles, the behavior of powder can be analyzed.
  • the behavior referred to here includes not only a change in position due to the movement of coarse-grained particles but also a change in state such as a temperature change.
  • the simulation method of the present embodiment may further include, for example, an initial setting step.
  • the positions of the particles constituting the powder to be analyzed can be initialized, and the analysis conditions, for example, the temperature of the region where the powder is arranged can be set if necessary.
  • the initial conditions are set in advance in the program used when analyzing the behavior of the coarse-grained particles in the coarse-grained particle behavior analysis step, or when the initial conditions are acquired by the first parameter acquisition step, the initial conditions are set.
  • the setting process does not have to be carried out.
  • the simulation method of the present embodiment may further include, for example, an output process.
  • the simulation result obtained in the coarse-grained particle behavior analysis step (S3) can be output to the output unit.
  • the content of the simulation result to be output is not particularly limited, but for example, the positions of the coarse-grained particles can be output and displayed as an image in time series in the output unit. Further, for example, the time-series change of the temperature distribution of the powder can be output as an image and displayed on the output unit.
  • the amount of calculation can be suppressed by forming a particle group composed of a plurality of particles into one coarse-grained particle. Therefore, it is possible to suppress the amount of calculation and efficiently perform the calculation even for the behavior of the powder on a large scale such as a plant used in a factory.
  • the simulation method of the present embodiment may further include a powder supply step, a reaction step, and the like in order to carry out various manufacturing steps using powder using the simulation results.
  • powder can be supplied from the powder supply device to the reaction furnace based on the result of the behavior of the powder in the reaction furnace obtained by the simulation. At this time, as the powder, a powder having predetermined physical properties selected based on the simulation result can be supplied.
  • the powder supplied to the reaction furnace in the powder supply process can be heated.
  • the powder can be heated under predetermined heating conditions based on the simulation result.
  • the powder having the desired physical properties is supplied from the powder supply device to the reaction furnace, and further determined in the reaction furnace.
  • the reaction ratio of the powder can be increased by heating under the heating conditions of. Further, since the heating conditions and the like can be optimized, the amount of energy used in the reaction can be optimized and the productivity can be improved.
  • the simulation method of the present embodiment implements the powder supply step or the like
  • the simulation method can also be referred to as a powder processing method or the like.
  • the program of this embodiment is related to a program for analyzing the behavior of a powder containing a plurality of particles, and the computer can function as each of the following parts.
  • the first parameter acquisition unit that acquires the first parameter including the parameters related to the powder.
  • a second parameter calculation unit that calculates the second parameter, which is a parameter for coarse-grained particles, when a particle group composed of a plurality of particles is coarse-grained and made into one coarse-grained particle.
  • Coarse-grained particle behavior analysis unit that analyzes the behavior of coarse-grained particles based on the first and second parameters.
  • the second parameter calculation unit can calculate the second parameter by using the solution of the characteristic equation using the relationship between the elastic energy of the particle group and the elastic energy of the coarse-grained particles.
  • the second parameter can include the thermal conductivity of the coarse-grained particles, and in this case, the second parameter calculation unit can calculate the thermal conductivity by using the solution of the above-mentioned characteristic equation.
  • the program of this embodiment can be stored in various storage media of a main storage device such as RAM or ROM of the simulation device described above or an auxiliary storage device, for example. Then, by reading and executing the program by the CPU, data can be read and written in the RAM or the like, and the input / output interface and the display device can be operated and executed. Therefore, the matters already explained in the simulation apparatus will be omitted.
  • the program of the present embodiment described above may be provided by storing it on a computer connected to a network such as the Internet and downloading it via the network. Further, the program of the present embodiment may be configured to be provided and distributed via a network such as the Internet.
  • the program of this embodiment may be distributed in a state of being stored in an optical disk such as a CD-ROM or a recording medium such as a semiconductor memory.
  • the amount of calculation can be suppressed by forming a particle group composed of a plurality of particles into one coarse-grained particle. Therefore, it is possible to suppress the amount of calculation and efficiently perform the calculation even for the behavior of the powder on a large scale such as a plant used in a factory.
  • the hat t in the following formula indicates the unit vector of the tangential overlap.
  • each amount of overlap between the particles and the wall surface is expressed by the following equation (29).
  • the subscript t in the following formula means a tangential component.
  • the amount of overlap between particles is expressed by the following formula (30).
  • the amount of overlap of the coarse-grained particles 21 is expressed by the following formula (31).
  • the vector ⁇ cw (indicated by an arrow above ⁇ cw in the following equation) is the rotation vector of the coarse-grained particles 21.
  • the degree of freedom in the rotation direction remains with respect to the movement of the center of gravity. Even if the positions of the centers of gravity match, the equation (31) and the equation (32) may not match. Therefore, it is necessary to convert the rotational motion of the particle group 11 and the rotational motion of the coarse-grained particles 21. Therefore, regarding rotation, it is assumed that the angular momentum and the rotational energy match before and after coarse graining. In this case, the angular momentum holds the following equation.
  • the first term and the second term on the right side mean the spinning motion component (Spin) and the revolution motion component (Orbit) of the particle group before coarse graining, respectively. ..
  • the vectors ⁇ cw , vector ⁇ s , and vector ⁇ o are of the coarse-grained particles. It means the angular velocity, the angular velocity of the rotation motion of the particle group before coarse graining, and the angular velocity of the revolving motion of the particle group before coarse graining.
  • I cw, I s, I o denotes the moment of inertia of the coarse-grained particles, rotation motion of the inertia moment of the coarse-grained front group of the particles, the coarse-grained inertia moment of the revolution of the front group of the particle.
  • the amount of overlap in the tangential direction of the coarse-grained particles can be defined as follows.
  • equations (33), (34), and equations (35) are used for simplification of calculation, but several types can be considered depending on the shape of the coarse-grained particles. For example, when trying to calculate the moment of inertia in consideration of the shape factor of the cubic coarse-grained particles, the following equations (38) to (41) are obtained.
  • Equation of motion in the rotational direction The torque in the tangential direction can be calculated by calculating the force in the tangential direction from the overlap amount in the tangential direction explained so far. Rolling friction between each particle and between a particle and a wall is said to generate a torque proportional to the product of the normal force applied to the particles and the contact radius. Since each rolling friction is generated in each particle, the rolling friction resistance of the entire coarse-grained particles can be expressed by the following equation (44).
  • a rolling frictional resistance of the entire vector T Tot_fric is coarse-grained particles
  • the vector T w a rolling frictional resistance between the Fric wall surface 12 or the external particles of the particle group 11
  • the vector T p, Fric Represents the rolling friction resistance between the internal particles of the particle group 11.
  • Each vector is indicated by an arrow above the character in the following formula.
  • the vectors T p and fric are doubled as shown in the equation (44).
  • the vectors T w, fric and the vectors T p, fric are represented by the equations (45) and (46).
  • the hat ⁇ in the equation means a unit vector in the rotation direction.
  • the contact radii r w and r p in the equations (45) and (46) can be geometrically calculated from the overlap amount, and are calculated based on the definition in the case of the above-mentioned tangential equation of motion. The accuracy can be improved by using the overlap amount.
  • ⁇ w means the rolling friction coefficient with the wall surface 12 of the particle group 11 or the external particles
  • ⁇ p means the rolling friction coefficient of the internal particles of the particle group 11.
  • the second parameter calculation unit uses the solution of the characteristic equation using the relationship between the elastic energy of the particle group and the elastic energy of the coarse-grained particles. Two parameters can be calculated. However, when the equation of motion in the tangential direction is used, the overlap amount can be calculated on the assumption that the angular momentum and the rotational energy match before and after coarse graining. Further, when the equation of motion in the rotation direction is used, the overlap amount calculated when the equation of motion in the tangential direction is used can be used when calculating the torque.
  • the simulation method of this embodiment the solution of the characteristic equation derived by using the relationship between the elastic energy of the grain group before coarse graining and the elastic energy of the coarse grained particles is used.
  • the second parameter can be calculated.
  • the overlap amount can be calculated on the assumption that the angular momentum and the rotational energy match before and after coarse graining.
  • the equation of motion in the rotation direction is used, the overlap amount calculated when the equation of motion in the tangential direction is used can be used when calculating the torque.
  • the second parameter calculation unit the solution of the characteristic equation derived by using the relationship between the elastic energy of the grain group before coarse graining and the elastic energy of the coarse grained particles is used.
  • the second parameter can be calculated.
  • the overlap amount can be calculated on the assumption that the angular momentum and the rotational energy match before and after coarse graining.
  • the overlap amount calculated when the equation of motion in the tangential direction is used can be used when calculating the torque.
  • [Comparative Example 1-2] The container was filled in a rectangular parallelepiped container under the same conditions as in Comparative Example 1-1 except that the coarse-grained magnification ⁇ was set to 2 and the particle group composed of 8 particles was made into 1 coarse-grained particle.
  • the temperature change of the powder layer was analyzed. Since the particles were coarse-grained, the particle size of the particles was changed as shown in Table 1, but the analysis was performed using the same parameters as in Comparative Example 1-1 except for the particle size.
  • FIG. 6 shows the change in the average temperature of the powder layer in the container obtained by the analysis.
  • the coarse-graining magnification ⁇ was set to 2, and the temperature change of the powder layer filled in the rectangular parallelepiped container was analyzed using the above-mentioned simulation device.
  • the first parameter acquisition unit 421 acquired the first parameter including the parameter related to the powder to be analyzed (first parameter acquisition step: S1).
  • the second parameter calculation unit 422 calculated the second parameter, which is a parameter for the coarse-grained particles (second parameter calculation step: S2).
  • Equation (13), Equation (19) Equation (26) using the parameter K r is a solution of the characteristic equation shown in, by the equation (27), the thermal conductivity of about coarse-grained particles Calculated.
  • the coefficient of restitution, the coefficient of friction, and the rolling friction coefficient were calculated using the characteristic equations. As described above, these coefficients are adjustable depending on the model applied to the calculation, and can be calculated and converted using the characteristic equations as described above.
  • the calculated parameters Kr are shown in Table 2.
  • Table 1 shows the first parameter and the second parameter obtained and calculated as described above. Then, based on the parameters shown in Table 1, the behavior of the coarse-grained particles, specifically, the temperature change was analyzed by the coarse-grained particle behavior analysis unit 423 (coarse-grained particle behavior analysis step: S3). The results are shown in FIG.
  • Example 1-1 are almost the same as the results of Comparative Example 1-1 without coarse graining, and the parameters of the coarse-grained particles are appropriately selected. , I was able to confirm that the analysis was performed.
  • Example 1-1 Since the coarse graining is performed in Example 1-1, the number of particles is 1/8 times that of Comparative Example 1-1, and the amount of calculation is suppressed as compared with Comparative Example 1-1. is made of.
  • Comparative Example 2 [Comparative Example 2-1] The behavior of the powder in the kiln, which is a rotating body, was analyzed by the discrete element method calculation using the parameters shown in Table 3. The coefficient of restitution was 0.75, the coefficient of friction was 0.3, and the coefficient of rolling friction was 0.5. The movement of the powder in the kiln obtained by analysis is shown in FIGS. 8A to 8D, the temperature distribution of the powder in the kiln obtained by analysis is shown in FIGS. 11A to 11C, and the powder (particles) in the kiln obtained by analysis is shown in FIGS. 11A to 11C. ) Are shown in FIG. 13, respectively.
  • FIG. 8A shows the state before the rotation of the kiln 70 is started, and the first powder group 71 and the second powder group 72 divided into groups are contained in half.
  • 8B, 8C, and 8D show the states when 2 seconds have passed, 4 seconds have passed, and 6 seconds have passed since the start of rotation of the kiln, respectively, and the first powder group 71 and the second powder group have passed. It shows a state in which 72 is mixed.
  • FIG. 11A shows the state before the rotation of the kiln is started, and it can be confirmed that the temperature is within the first temperature range 101 and is uniform.
  • 11B and 11C show the states 3 seconds and 6 seconds after the start of rotation of the kiln, respectively. Since heating is performed from the outer wall side of the kiln 70, the powder is heated in order from the center side. It can be confirmed that the first temperature range 101, the second temperature range 102, and the third temperature range 103 are distributed. The temperature increases in the order of the first temperature range 101, the second temperature range 102, and the third temperature range 103.
  • FIG. 9A shows the state before the rotation of the kiln 70 is started, and the first powder group 71 and the second powder group 72 divided into groups are contained in half.
  • 9B, 9C, and 9D show the states when 2 seconds have passed, 4 seconds have passed, and 6 seconds have passed since the start of rotation of the kiln, respectively, and the first powder group 71 and the second powder group have passed. It shows a state in which 72 is mixed.
  • FIG. 12A shows the state before the rotation of the kiln is started, and it can be confirmed that the temperature is within the first temperature range 101 and is uniform.
  • 12B and 12C show the states 3 seconds and 6 seconds after the start of rotation of the kiln, respectively. Since heating is performed from the outer wall side of the kiln 70, the powder is heated in order from the center side. It can be confirmed that the first temperature range 101, the second temperature range 102, and the third temperature range 103 are distributed. The temperature increases in the order of the first temperature range 101, the second temperature range 102, and the third temperature range 103.
  • Example 2-1 The coarse-graining magnification ⁇ was set to 4, and the behavior of the powder in the kiln, which is a rotating body, was analyzed using the simulation device described above.
  • the first parameter acquisition unit 421 acquired the first parameter including the parameter related to the powder to be analyzed (first parameter acquisition step: S1).
  • the second parameter calculation unit 422 calculated the second parameter, which is a parameter for the coarse-grained particles (second parameter calculation step: S2).
  • Equation (13), Equation (19) Equation (26) using the parameter K r is a solution of the characteristic equation shown in, by the equation (27), the thermal conductivity of about coarse-grained particles Calculated.
  • the coefficient of restitution, the coefficient of friction, and the rolling friction coefficient were calculated using the characteristic equations.
  • the calculated parameters Kr are shown in Table 4.
  • Table 3 shows the first parameter and the second parameter acquired and calculated as described above. Then, based on the parameters shown in Table 3, the behavior of the coarse-grained particles, specifically, the movement in the kiln and the temperature change were analyzed by the coarse-grained particle behavior analysis unit 423 (coarse-grained particle behavior analysis). Step: S3).
  • the movement of the powder in the kiln obtained by analysis is shown in FIGS. 7A to 7D
  • the temperature distribution of the powder in the kiln obtained by analysis is shown in FIGS. 10A to 10C
  • the average of the powder in the kiln obtained by analysis is shown in FIGS. 10A to 10C.
  • the temperatures are shown in FIG. 13, respectively.
  • FIG. 7A shows the state before the rotation of the kiln 70 is started, and the first powder group 71 and the second powder group 72 divided into groups are contained in half.
  • 7B, 7C, and 7D show the states when 2 seconds have passed, 4 seconds have passed, and 6 seconds have passed since the start of rotation of the kiln, respectively, and the first powder group 71 and the second powder group have passed. It shows a state in which 72 is mixed.
  • FIG. 10A shows the state before the rotation of the kiln is started, and it can be confirmed that the temperature is within the first temperature range 101 and is uniform.
  • FIGS. 10B and 10C show the states 3 seconds and 6 seconds after the start of rotation of the kiln, respectively. Since heating is performed from the outer wall side of the kiln 70, the powder is heated in order from the center side. It can be confirmed that the first temperature range 101, the second temperature range 102, and the third temperature range 103 are distributed. The temperature increases in the order of the first temperature range 101, the second temperature range 102, and the third temperature range 103.
  • Example 2-1 are not coarse-grained. It was confirmed that the results of Comparative Example 2-1 were almost the same, and that the parameters of the coarse-grained particles were appropriately selected and set, and the analysis was performed.
  • Example 2-1 and Comparative Example 2-1 showed the same tendency such as having an inflection point.
  • Example 2-1 Since the coarse graining is performed in Example 2-1 the number of particles is 1/64 times that of Comparative Example 2-1 and the amount of calculation is calculated as compared with Comparative Example 2-1. Can be suppressed.

Landscapes

  • Engineering & Computer Science (AREA)
  • Chemical & Material Sciences (AREA)
  • Theoretical Computer Science (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Computing Systems (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Physics & Mathematics (AREA)
  • General Health & Medical Sciences (AREA)
  • Health & Medical Sciences (AREA)
  • Crystallography & Structural Chemistry (AREA)
  • General Physics & Mathematics (AREA)
  • Biochemistry (AREA)
  • Pathology (AREA)
  • Immunology (AREA)
  • Dispersion Chemistry (AREA)
  • Analytical Chemistry (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Medicinal Chemistry (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

複数の粒子を含む粉体の挙動を解析するためのシミュレーション装置であって、　前記粉体に関連するパラメータを含む第１パラメータを取得する第１パラメータ取得部と、　複数個の前記粒子から構成される粒子群を粗視化し、１個の粗視化粒子とした場合の、前記粗視化粒子についてのパラメータである第２パラメータを算出する第２パラメータ算出部と、　前記第１パラメータ、および前記第２パラメータに基づいて、前記粗視化粒子の挙動を解析する粗視化粒子挙動解析部と、を有し、　前記第２パラメータ算出部は、前記粒子群の弾性エネルギーと、前記粗視化粒子の弾性エネルギーとの関係を用いた特性方程式の解を用いて、前記第２パラメータを算出するシミュレーション装置を提供する。

Description

シミュレーション装置、シミュレーション方法、プログラム
　本発明は、シミュレーション装置、シミュレーション方法、プログラム
に関する。
　特許文献１には、シミュレーション条件の入力を行う入力装置と、
　シミュレーション結果を出力する出力装置と、
　前記入力装置から入力されたシミュレーション条件に基づいて、大きさが異なる複数の
粒子を含む粉粒体の挙動を解析する処理装置とを有し、
　前記処理装置は、
　前記入力装置から入力されたシミュレーション対象の粉粒体の粒径分布を規定するパラメータの値、及び粒子を粗視化する基準となる粗視化係数の値に基づいて、粗視化された粉粒体の挙動をシミュレーションにより求め、
　シミュレーションによって求められた粒子の挙動と、入力された粗視化係数の値とを関連付けて前記出力装置に出力するシミュレーション装置等に開示されているように、粉粒体についてのシミュレーション装置が開示されている。
日本国特開２０２０－５７１３５号公報
　工場等におけるプロセスの改善や、製造工程を検討する際の試験工数の削減等を目的として、複数の粒子を含む粉体（粉粒体）の挙動の解析を離散要素法（ＤＥＭ：Ｄｉｓｃｒｅｔｅ　Ｅｌｅｍｅｎｔ　Ｍｅｔｈｏｄ）計算等により行うことが従来からなされてきた。
　離散要素法計算は、個々の粒子について運動方程式を解くことで、粉体全体の運動を記述するシミュレーション技術である。
　しかしながら、離散要素法計算では、取り扱う粒子数が多くなるほど計算負荷が高くなる。このため、例えば工場で用いるプラント等のような大きな規模での粉体の挙動について解析を行う場合、計算量が膨大になり、現実的には計算を行うことが困難となる。
　そこで、複数個の粒子から構成される粒子群を１個の粗視化粒子とする粗視化の手法を用いたシミュレーション装置について従来から検討されていた（特許文献１を参照）。粗視化の手法を用いたシミュレーション装置では、正確な解析結果を得るため、計算で用いる粗視化粒子に関するパラメータを適切に選択することが求められる。このため、粗視化粒子に関するパラメータを新たな手法により選択、設定し、複数の粒子を含む粉体の挙動を解析できる新たなシミュレーション装置が求められていた。
　上記従来技術が有する問題に鑑み、本発明の一側面では、複数の粒子を含む粉体の挙動を解析できる新たなシミュレーション装置を提供することを目的とする。
　上記課題を解決するため本発明の一態様によれば、
　複数の粒子を含む粉体の挙動を解析するためのシミュレーション装置であって、
　前記粉体に関連するパラメータを含む第１パラメータを取得する第１パラメータ取得部と、
　複数個の前記粒子から構成される粒子群を粗視化し、１個の粗視化粒子とした場合の、前記粗視化粒子についてのパラメータである第２パラメータを算出する第２パラメータ算出部と、
　前記第１パラメータ、および前記第２パラメータに基づいて、前記粗視化粒子の挙動を解析する粗視化粒子挙動解析部と、を有し、
　前記第２パラメータ算出部は、前記粒子群の弾性エネルギーと、前記粗視化粒子の弾性エネルギーとの関係を用いた特性方程式の解を用いて、前記第２パラメータを算出するシミュレーション装置を提供する。
　本発明の一態様によれば、複数の粒子を含む粉体の挙動を解析できる新たなシミュレーション装置を提供することができる。
複数個の粒子から構成される粒子群と、粒子群と壁面との衝突時の説明図である。 複数個の粒子から構成される粒子群と、粒子群と壁面との衝突時の説明図である。 複数個の粒子から構成される粒子群を粗視化した粗視化粒子と、粗視化粒子と壁面との衝突時の説明図である。 複数個の粒子から構成される粒子群を粗視化した粗視化粒子と、粗視化粒子と壁面との衝突時の説明図である。 本発明の実施形態に係るシミュレーション装置のハードウェア構成図である。 本発明の実施形態に係るシミュレーション装置の機能を示すブロック図である。 本発明の実施形態に係るシミュレーション方法を説明するフローチャートである。 実験例１における粉体の平均温度の変化を示すグラフである。 実施例２－１におけるキルン内での粉体の混合の状態を示す図である。 実施例２－１におけるキルン内での粉体の混合の状態を示す図である。 実施例２－１におけるキルン内での粉体の混合の状態を示す図である。 実施例２－１におけるキルン内での粉体の混合の状態を示す図である。 比較例２－１におけるキルン内での粉体の混合の状態を示す図である。 比較例２－１におけるキルン内での粉体の混合の状態を示す図である。 比較例２－１におけるキルン内での粉体の混合の状態を示す図である。 比較例２－１におけるキルン内での粉体の混合の状態を示す図である。 比較例２－２におけるキルン内での粉体の混合の状態を示す図である。 比較例２－２におけるキルン内での粉体の混合の状態を示す図である。 比較例２－２におけるキルン内での粉体の混合の状態を示す図である。 比較例２－２におけるキルン内での粉体の混合の状態を示す図である。 実施例２－１におけるキルン内の粉体の温度分布を示す図である。 実施例２－１におけるキルン内の粉体の温度分布を示す図である。 実施例２－１におけるキルン内の粉体の温度分布を示す図である。 比較例２－１におけるキルン内の粉体の温度分布を示す図である。 比較例２－１におけるキルン内の粉体の温度分布を示す図である。 比較例２－１におけるキルン内の粉体の温度分布を示す図である。 比較例２－２におけるキルン内の粉体の温度分布を示す図である。 比較例２－２におけるキルン内の粉体の温度分布を示す図である。 比較例２－２におけるキルン内の粉体の温度分布を示す図である。 実験例２における粉体の平均温度の変化を示すグラフである。
　本開示の一実施形態（以下「本実施形態」と記す）に係るシミュレーション装置、シミュレーション方法、プログラムの具体例を、以下に図面を参照しながら説明する。なお、本発明はこれらの例示に限定されるものではなく、特許請求の範囲によって示され、特許請求の範囲と均等の意味および範囲内での全ての変更が含まれることが意図される。
１．第１の実施形態
［シミュレーション装置］
（１）粒子の粗視化、および粗視化粒子の粒子挙動の計算に用いるパラメータについて
（１－１）粒子の粗視化について
　本実施形態のシミュレーション装置の詳細について説明する前に、本実施形態のシミュレーション装置で用いることができる、複数個の粒子から構成される粒子群の粗視化、および粗視化した粒子である粗視化粒子に関連するパラメータの算出方法について以下に説明する。
　既述のように、離散要素法計算では、取り扱う粒子数が多くなるほど計算負荷が高くなる。このため、例えば工場で用いるプラント等のような大きな規模での粉体の挙動について解析を行う場合、計算量が膨大になり、現実的には計算を行うことが困難となる。
　そこで、多量の粒子を含む粉体について挙動を解析する場合、計算量を抑制するため、例えば図１Ａに示した複数個の粒子から構成される粒子群１１を、図２Ａに示すような１個の大きな粒子である粗視化粒子２１として取り扱う粗視化の技術が必要となる。
　ただし、粗視化前の個別の粒子と、粗視化粒子とでは比表面積等が異なるため、計算に要する一部のパラメータは変化することになる。このため、粗視化粒子のパラメータを適切に決定する必要がある。
（１－２）粗視化粒子の粒子挙動の計算に用いるパラメータについて
　そこで、本発明の発明者は、粗視化粒子に関するパラメータを決定する方法について検討を行った。計算にあたって、図１Ａに示した粗視化する前の複数個の粒子から構成される粒子群１１が壁面１２に衝突する場合と、図２Ａに示した粗視化粒子２１が壁面１２に衝突する場合とをモデルに用いた。以下の説明では壁面と粒子とが衝突する場合を例に粗視化粒子に関するパラメータを決定する方法を記載するが、粒子同士が衝突する場合でも同様の議論となるため説明を省略する。
　図１Ａに示すように、複数個の粒子から構成される粒子群１１が、立方体状に、縦方向、横方向、高さ方向に２個ずつ、合計で２３個並んでいるとする。後述するように、係る８個の粒子をまとめて１個の粗視化粒子とする場合、一辺の方向に並べた粒子の数、すなわち２を粗視化倍率とする。
　図１Ａに示した複数個の粒子１１Ａ、１１Ｂを含む粒子群１１が壁面１２に衝突する際に、粒子群１１のうち壁面１２側に位置する粒子１１Ａが壁面あるいは外部粒子から受ける力を図１Ｂに示すようにＦｗとする。また、図１Ｂに示すように粒子１１Ａの壁面１２あるいは外部粒子とのオーバーラップ量をδｗ、粒子１１Ｂの隣接する粒子１１Ａとのオーバーラップ量をδｐとする。なお、図１Ｂは、粒子群１１が壁面１２に衝突する際の様子を側面側から見た図になる。
　この場合、粒子群１１に加わる力の大きさは、以下の式（１）で表すことができる。
　なお、式（１）におけるαは粗視化倍率であり、既述のように粒子群１１を１個の粗視化粒子とした場合に一辺の方向に並べた粒子の数を意味する。図１Ａに示した粒子群１１を図２Ａに示した１個の粗視化粒子２１とする場合、α＝２となる。
　また、ｍは各粒子１１Ａ、１１Ｂの質量を、ａＧは粒子群１１の重心の加速度を、ηｗは壁面１２あるいは外部粒子と粒子１１Ａとの反発係数から算出される粘性係数をそれぞれ意味している。作用反作用の法則により粒子間の接触力が打ち消されるため、壁面１２とは直接接しない粒子１１Ｂが、隣接する粒子１１Ａから受ける力Ｆｐは式（１）には表れないことになる。
Figure JPOXMLDOC01-appb-M000001
 　なお、上述のηｗ等を反発係数から算出する際に用いる、反発係数ｅと粘性係数ηとの関係は、以下の式（Ａ）で表すことができる。式（Ａ）中のｍ＊は換算質量を、Ｋはバネ係数を意味している。
Figure JPOXMLDOC01-appb-M000002
 　次に、図２Ａに示すように、図１Ａに示した８個の粒子からなる粒子群１１を１つの粗視化粒子２１としたとする。この場合において、粗視化粒子２１が、壁面１２に衝突した場合に、粗視化粒子２１が受ける力は以下の式（２）で表すことができる。
　式（２）中のＦｃｗは、図２Ｂに示すように粗視化粒子２１が壁面１２あるいは外部粒子から受ける力を、δｃｗは粗視化粒子２１の壁面１２あるいは外部粒子とのオーバーラップ量を、ηｃｗは粗視化粒子２１と壁面１２あるいは外部粒子との反発係数から算出される粘性係数をそれぞれ意味している。なお、図２Ｂは、粗視化粒子２１が壁面１２に衝突する際の様子を側面側から見た図になる。
Figure JPOXMLDOC01-appb-M000003
 　既述のように、粗視化は、離散要素法計算において、計算量を抑制するために行うものである。このため、粗視化粒子２１についての計算結果と、粗視化粒子２１とする前の粒子群１１についての計算結果とは一致することになる。
　従って、粒子群１１について計算を行った既述の式（１）と、粒子群を粗視化した粗視化粒子２１について計算を行った既述の式（２）とから、対応するパラメータが一致することを示す以下の式（３）、式（４）が導出される。
Figure JPOXMLDOC01-appb-M000004
Figure JPOXMLDOC01-appb-M000005
 　また、Ｈｅｒｔｚ－Ｍｉｎｄｌｉｎの接触モデルにより、粒子のオーバーラップ量δｗ、δｐ、δｃｗを用いて、各粒子に加わる力を以下の式（５）～式（７）のように表すことができる。式（５）～式（７）におけるＫｗは粒子１１Ａと壁面１２あるいは外部粒子とのバネ係数を、Ｋｐは粒子群１１の内部粒子のバネ係数を、Ｋｃｗは粗視化粒子２１と壁面１２あるいは外部粒子とのバネ係数をそれぞれ意味している。
Figure JPOXMLDOC01-appb-M000006
Figure JPOXMLDOC01-appb-M000007
Figure JPOXMLDOC01-appb-M000008
 　そうすると、式（３）、式（５）、式（７）から以下の式（８）の関係が導かれる。
Figure JPOXMLDOC01-appb-M000009
 　ここで、以下の式（９）のようにＫｒを定義すると、上記式（８）は以下の式（１０）のように表すことができる。
Figure JPOXMLDOC01-appb-M000010
Figure JPOXMLDOC01-appb-M000011
 　なお、粗視化前の粒子群１１と、粗視化粒子２１とで重心が一致しているとすると、以下の式（１１）の関係を満たすことになる。
Figure JPOXMLDOC01-appb-M000012
 　従って、Ｋｒを適切に設定することで、粗視化粒子２１の壁面１２あるいは外部粒子とのオーバーラップ量δＣＷから、粗視化前の粒子群１１を構成する粒子のオーバーラップ量を算出できることもわかる。
　そして、Ｋｒは、衝突中の粗視化前の粒子群１１の弾性エネルギーと、粗視化粒子２１の弾性エネルギーとの関係を用いた特性方程式により算出できる。具体的には例えば、粗視化前の粒子群１１全体の弾性エネルギーと、粗視化粒子全体の弾性エネルギーが等しくなると仮定して、特性方程式を作成してＫｒを算出できる。
　壁面１２と衝突中の、粒子群１１の弾性エネルギーと、粗視化粒子の弾性エネルギーは、既述の粒子群１１を構成する粒子１１Ａ、１１Ｂに加わる力や、粗視化粒子に加わる力を表した式（５）～式（７）をオーバーラップ量の距離で積分することで算出できる。
　従って、粗視化前の粒子群１１全体の弾性エネルギーと、粗視化粒子全体の弾性エネルギーとを用いて、以下の式（１２）が得られる
Figure JPOXMLDOC01-appb-M000013
 　上記式（１２）は、既述の式（８）～式（１１）を用いて以下の式（１３）に変形できる。
Figure JPOXMLDOC01-appb-M000014
 　式（１３）は垂直方向のＫｒの特性方程式である。そして、式（９）の定義式から明らかなように、Ｋｒは粗視化粒子、および粗視化前の粒子群１１を構成する粒子のオーバーラップ量に関するパラメータであり、粗視化粒子の挙動を支配するパラメータである。このため、特性方程式によりＫｒを事前に求めておくことによって、粗視化後粒子のオーバーラップ量から粗視化前粒子群のオーバーラップ量を算出するなど、粗視化粒子についてのパラメータを算出したり、粗視化粒子の挙動を計算したりすることが可能になる。
　ここまで、壁面１２に対する垂直方向の運動方程式を用いて説明したが、接線方向の運動方程式や、回転の運動方程式についても同様のことが言える。
　具体的には、接線方向の運動方程式は、式（１４）で表すことができる。
　この場合も、式（１５）に示したように、Ｋｒを設定すると、δｗ、δｐは式（１６）、式（１７）のように示すことができ、粗視化前の粒子群の弾性エネルギーと、粗視化粒子の弾性エネルギーが等しいとすると、式（１８）が得られる。式（１８）を変形することで、接線方向の特性方程式である式（１９）が得られる。ただし、接線方向の接触モデルには線形バネモデルを使用した。このように、接触モデルによって弾性エネルギーの計算式が異なるが、必要に応じて特性方程式を変更する等により、適切に弾性エネルギーを算出できる。
Figure JPOXMLDOC01-appb-M000015
Figure JPOXMLDOC01-appb-M000016
Figure JPOXMLDOC01-appb-M000017
Figure JPOXMLDOC01-appb-M000018
Figure JPOXMLDOC01-appb-M000019
Figure JPOXMLDOC01-appb-M000020
 　回転の運動方程式については第２の実施形態で詳述する。
（１－３）粗視化粒子の伝熱に関するパラメータについて
　既述のＫｒを用いて、粗視化粒子の伝熱に関するパラメータである熱伝導率を求めることもできる。
　粒子の伝熱は熱伝導率を用いて、以下の式（２０）、式（２１）のようにあらわされる。
　なお、式（２０）、式（２１）において、Ｑは熱流量を、ｈは熱流量係数を、ΔＴは壁面あるいは粒子同士の温度差を、ｋｗは壁面１２あるいは外部粒子と粒子１１Ａとの熱伝導率を、ｋｐは粒子群１１内部の熱伝導率を、ａは粒子群１１の接触半径をそれぞれ意味する。
Figure JPOXMLDOC01-appb-M000021
Figure JPOXMLDOC01-appb-M000022
 　粒子と壁面や、他の粒子との接触半径は、粒径によって変化するため、粒径が大きくなると粒子に与えられる熱流量が変化し、粒子の温度変化に影響する。
　接触半径は、各粒子のオーバーラップ量に依存する。そして、既述のように粗視化前後でのオーバーラップ量は既述の垂直方向の特性方程式の解に関係づけられている。このことから、垂直方向の特性方程式の解Ｋｒを用いて伝熱方程式を書き換えると以下の式（２２）、式（２３）のようになる。
　式（２３）におけるａ´は粗視化粒子の接触半径を、ｒは粒子群１１を構成する粒子の半径をそれぞれ意味する。
Figure JPOXMLDOC01-appb-M000023
Figure JPOXMLDOC01-appb-M000024
 　そして、粗視化粒子の伝熱方程式は以下の式（２４）、式（２５）で表すことができる。
　式（２４）、式（２５）におけるＱｃは粗視化粒子の熱流量を、ｈｃは粗視化粒子の熱流量係数を、ΔＴｃは粗視化粒子２１と壁面１２あるいは外部粒子との温度差を、ｋｗ´は壁面１２と粗視化粒子２１との熱伝導率を、ｋｐ´は外部粒子と粗視化粒子との熱伝導率をそれぞれ意味する。
Figure JPOXMLDOC01-appb-M000025
Figure JPOXMLDOC01-appb-M000026
 　熱伝導率に注目すると、粗視化前後での熱伝導率を、Ｋｒを用いて以下の式（２６）、式（２７）のように変換すると、粗視化前後で等価な熱伝導方程式を得ることができる。ただし、ΔＴｃ＝αΔＴとした。
Figure JPOXMLDOC01-appb-M000027
Figure JPOXMLDOC01-appb-M000028
 　つまり粗視化する際に、熱伝導率をＫｒ １／３α１／２倍することで、粗視化前の粒子群に与えられる熱流量と、粗視化後の粒子に与えられる熱流量を等しくすることができ、その結果、温度の時間変化を粗視化前後で一致させることができる。
　なお、ここでは粗視化前後での熱伝導率を例に説明したが、熱伝導率以外のパラメータについても同様にして、既述の特性方程式の解Ｋｒを用いて粗視化後のパラメータを算出できる。例えば、特性方程式を用いて反発係数、摩擦係数、転がり摩擦係数を算出することもできる。なお、計算に適用するモデルに応じて、これらの係数は調整可能であり、上述のように特性方程式を用いて算出、換算できる。
（２）シミュレーション装置
　本実施形態のシミュレーション装置は、複数の粒子を含む粉体の挙動を解析するためのシミュレーション装置であり、以下の部材を有することができる。　
　粉体に関連するパラメータを含む第１パラメータを取得する第１パラメータ取得部。　
　複数個の粒子から構成される粒子群を粗視化し、１個の粗視化粒子とした場合の、粗視化粒子についてのパラメータである第２パラメータを算出する第２パラメータ算出部。　
　第１パラメータ、および第２パラメータに基づいて、粗視化粒子の挙動を解析する粗視化粒子挙動解析部。
　第２パラメータ算出部は、粒子群の弾性エネルギーと、粗視化粒子の弾性エネルギーとの関係を用いた特性方程式の解を用いて、第２パラメータを算出する。
　図３に示したハードウェア構成図に示すように、本実施形態のシミュレーション装置３０は、例えば、情報処理装置（コンピュータ）で構成され、物理的には、演算処理部であるＣＰＵ（Ｃｅｎｔｒａｌ　Ｐｒｏｃｅｓｓｉｎｇ　Ｕｎｉｔ：プロセッサ）３１と、主記憶装置であるＲＡＭ（Ｒａｎｄｏｍ　Ａｃｃｅｓｓ　Ｍｅｍｏｒｙ）３２やＲＯＭ（Ｒｅａｄ　Ｏｎｌｙ　Ｍｅｍｏｒｙ）３３と、補助記憶装置３４と、入出力インタフェース３５と、出力装置である表示装置３６等を含むコンピュータシステムとして構成することができる。これらは、バス３７で相互に接続されている。なお、補助記憶装置３４や表示装置３６は、外部に設けられていてもよい。
　ＣＰＵ３１は、シミュレーション装置３０の全体の動作を制御し、各種の情報処理を行う。ＣＰＵ３１は、ＲＯＭ３３または補助記憶装置３４に格納された、例えば後述するシミュレーション方法や、プログラム（シュミレーションプログラム）を実行して、粗視化粒子についてのパラメータである第２パラメータの算出や、粗視化粒子の挙動の解析を行うことができる。
　ＲＡＭ３２は、ＣＰＵ３１のワークエリアとして用いられ、主要な制御パラメータや情報を記憶する不揮発ＲＡＭを含んでもよい。
　ＲＯＭ３３は、プログラム（シュミレーションプログラム）等を記憶することができる。
　補助記憶装置３４は、ＳＳＤ（Ｓｏｌｉｄ　Ｓｔａｔｅ　Ｄｒｉｖｅ）や、ＨＤＤ（Ｈａｒｄ　Ｄｉｓｋ　Ｄｒｉｖｅ）等の記憶装置であり、シミュレーション装置の動作に必要な各種のデータ、ファイル等を格納できる。
　入出力インタフェース３５は、タッチパネル、キーボード、表示画面、操作ボタン等のユーザインタフェースと、外部のデータ収録サーバ等からの情報を取り込み、他の電子機器に解析情報を出力する通信インタフェースとの双方を含む。
　表示装置３６は、モニタディスプレイ等である。表示装置３６では、解析画面が表示され、入出力インタフェース３５を介した入出力操作に応じて画面が更新される。
　図３に示したシミュレーション装置３０の各機能は、例えばＲＡＭ３２やＲＯＭ３３等の主記憶装置または補助記憶装置３４からプログラム（シミュレーションプログラム）等を読み込ませ、ＣＰＵ３１により実行することにより、ＲＡＭ３２等におけるデータの読み出しおよび書き込みを行うと共に、入出力インタフェース３５および表示装置３６を動作させることで実現できる。
　図４に、本実施形態のシミュレーション装置３０の機能ブロック図を示す。
　図４に示すように、シミュレーション装置３０は、受付部４１、処理装置４２、出力部４３を有することができる。これらの各部は、シミュレーション装置３０が有するＣＰＵ、記憶装置、各種インタフェース等を備えたパーソナルコンピュータ等の情報処理装置において、ＣＰＵが予め記憶されている例えば後述するシミュレーション方法や、プログラムを実行することでソフトウェアおよびハードウェアが協働して実現される。
　各部の構成について以下に説明する。
（Ａ）受付部
　受付部４１は、処理装置４２で実行される処理に関係するユーザーからのコマンドやデータの入力を受け付ける。受付部４１としてはユーザーが操作を行い、コマンド等を入力するキーボードやマウス、ネットワークを介して入力を行う通信装置、ＣＤ－ＲＯＭ、ＤＶＤ－ＲＯＭ等の各種記憶媒体から入力を行う読み取り装置などが挙げられる。
（Ｂ）処理装置
　処理装置４２は、第１パラメータ取得部４２１、第２パラメータ算出部４２２、粗視化粒子挙動解析部４２３を有することができる。なお、処理装置は、必要に応じてさらに任意の部材を有することができ、例えば初期設定部等を有することもできる。
（Ｂ－１）第１パラメータ取得部
　第１パラメータ取得部４２１では、例えば解析の対象となる粉体に関連するパラメータを含む第１パラメータを取得できる。第１パラメータは、粉体に関連するパラメータ以外に、解析に要する各種パラメータを含むこともできる。第１パラメータは解析（シミュレーション）の内容に応じて選択できるため、その具体的な種類は特に限定されない。第１パラメータとしては、離散要素法計算に必要となる各種パラメータが挙げられ、具体的には例えば粒子径、粒子数、ヤング率、計算のＴｉｍｅ　ｓｔｅｐ、ポワソン比、壁面との摩擦係数、粒子間の摩擦係数、転がり摩擦係数、密度等から選択された１種類以上が挙げられる。
　第１パラメータは、データベース等に収録されているデータであってもよく、予め実験を行い求めた実験値であってもよい。また、第１パラメータは、実験結果からシミュレーション等によりフィッティングして算出した計算値であってもよい。
（Ｂ－２）第２パラメータ算出部
　既述のように、本実施形態のシミュレーション装置３０では、計算量を抑制するため、粉体が有する複数個の粒子から構成される粒子群を、１個の粗視化粒子に粗視化し、粒子の数を少なくして計算を行うことができる。ただし、粗視化粒子は、粗視化前の粒子群を構成する個々の粒子とは質量等の各種パラメータが異なる。このため、粗視化粒子について、計算を行う際に必要となるパラメータの算出や、設定を行う必要がある。
　第２パラメータ算出部４２２では、「（１）粒子の粗視化、および粗視化粒子の粒子挙動の計算に用いるパラメータについて」の中で説明したように、粗視化前の粒子群の弾性エネルギーと、粗視化粒子の弾性エネルギーとの関係を用いて導出した特性方程式の解であるＫｒを用いて、第２パラメータを算出できる。具体的には例えば、粗視化前の粒子群１１全体の弾性エネルギーと、粗視化粒子全体の弾性エネルギーが等しくなると仮定して、既述の垂直方向のＫｒの特性方程式である式（１３）を導出し、式（１３）から垂直方向のＫｒを算出できる。そして、既述の式（１３）に示した特性方程式の解である垂直方向のＫｒを用いて、第２パラメータを算出できる。また、既述の式（１９）に示した接線方向のＫｒの特性方程式を用いて、接線方向のＫｒを算出でき、係る接線方向のＫｒを用いて第２パラメータを算出することもできる。
　既述のように、Ｋｒは、粗視化粒子の挙動を支配するパラメータであり、Ｋｒを用いることで、粗視化粒子の挙動に関連する各種パラメータを算出できる。
　後述する粗視化粒子挙動解析部で用いる第２パラメータの種類は解析の内容に応じて選択できるため、特に限定されない。例えば第２パラメータは、粗視化粒子の熱伝導率を含むこともできる。この場合、第２パラメータ算出部は、既述の特性方程式の解Ｋｒを用いて、熱伝導率を算出できる。
（Ｂ－３）粗視化粒子挙動解析部
　粗視化粒子挙動解析部４２３では、第１パラメータ取得部４２１により取得した第１パラメータ、および第２パラメータ算出部４２２で算出した第２パラメータを用いて、粗視化粒子の挙動を解析できる。具体的には離散要素法を用いて計算し、粗視化粒子の挙動を解析できる。粗視化粒子の挙動を解析することで、粉体の挙動を解析することができる。
　なお、ここでいう挙動とは、粗視化粒子の運動による位置の変化だけではなく、温度変化等の状態変化も含む。
（Ｂ－４）初期設定部
　図示しない初期設定部は、解析対象となる粉体を構成する粒子の位置を初期化するとともに、解析の条件、例えば必要に応じて粉体を配置する領域の温度等を設定できる。なお、例えば粗視化粒子挙動解析部４２３で粗視化粒子の挙動を解析する際に用いるプログラム等に予め初期条件が設定されている場合や、第１パラメータ取得部４２１により取得する場合には、初期設定部は設けなくてもよい。
（Ｃ）出力部
　出力部４３は、ディスプレイ等を有することができる。粗視化粒子挙動解析部４２３で得られたシミュレーション結果を出力部４３に出力できる。出力するシミュレーション結果の内容は特に限定されないが、例えば出力部４３に粗視化された粒子の位置を時系列で画像として出力し、表示することができる。また、例えば出力部４３に、粉体の温度分布の時系列変化を画像として出力し、表示することもできる。
　以上に説明した本実施形態のシミュレーション装置によれば、複数の粒子を含む粉体の挙動をシミュレーションすることができ、その用途等は特に限定されない。例えば、キルン等の回転体内での粉体の挙動のシミュレーションに好適に用いることができる。すなわち、本実施形態のシミュレーション装置によれば、回転体内での、粉体の挙動を解析することもできる。
　以上に説明した本実施形態のシミュレーション装置によれば、複数個の粒子からなる粒子群を１個の粗視化粒子とすることで、計算量を抑制できる。このため、工場で用いるプラント等のような大きな規模での粉体の挙動についても計算量を抑制し、効率的に計算を行うことができる。
　そして、粗視化粒子のパラメータを既述のパラメータＫｒを用いて計算するため、精度よく計算を行うことができる。
　なお、本実施形態のシミュレーション装置は、シミュレーション結果を用いて、粉体を用いた各種製造工程を実施するために、さらに粉体供給装置、反応炉、制御装置等を有することもできる。
　粉体供給装置としては、ホッパー等の粉体を貯蔵、排出できる装置が挙げられる。粉体供給装置は、ホッパー等から反応炉へと排出、供給する粉体の量を制御するため、フィーダーや、弁等の供給量調整装置をさらに有することもできる。シミュレーション結果に基づいて、所望の物性、例えば所望の平均粒径の粉体を供給できるように、異なる物性の粉体をそれぞれ収容した複数の粉体供給装置を有することが好ましい。
　反応炉としては、加熱炉等の各種反応炉が挙げられ、回転炉、例えばキルン等が挙げられる。粉体供給装置と、反応炉とは配管により接続しておくことができる。
　制御装置は、得られた反応炉内での粉体の挙動の結果に基づいて、粉体供給装置に対して、反応炉に所望の物性、例えば所望の平均粒径の粉体を供給するように制御できる。また、制御装置は、得られた反応炉内での粉体の挙動の結果に基づいて、反応炉の加熱条件を制御できる。加熱条件としては、例えば反応炉内の温度条件や、雰囲気条件、加熱時間等が挙げられる。なお、粉体供給装置や反応炉には、粉体の供給量や、温度等を検出するために各種センサを設けておき、任意のタイミングで測定した測定結果を制御装置に供給できるように構成しても良い。この場合、制御装置は得られた測定データに基づいて各装置の制御を行うこともできる。
　制御装置は、制御条件等のデータ処理を行い、粉体供給装置や、反応炉との通信が行えるように、ＣＰＵや、主記憶装置、補助記憶装置、入出力インタフェース等を有することができる。主記憶装置としてはＲＡＭや、ＲＯＭが、補助記憶装置としてはＳＳＤや、ＨＤＤ等が挙げられる。入出力インタフェースとしては、粉体供給装置や、反応炉との間で制御信号や、データをやり取りするための通信インタフェースが挙げられる。通信インタフェースの種類は特に限定されない。有線、無線のいずれの通信方法も用いることができ、例えば有線ＬＡＮ（Ｌｏｃａｌ　Ａｒｅａ　Ｎｅｔｗｏｒｋ）や、無線ＬＡＮ等が挙げられる。
　上述のように、シミュレーションにより得られた反応炉内での粉体の挙動の結果に基づいて、粉体供給装置から所望の物性の粉体を供給し、さらに反応炉内で所定の加熱条件で加熱することで、粉体の反応割合を高めることができる。また、加熱条件等を最適化できるため、反応の際に用いるエネルギー量を最適化し、生産性を向上できる。
　なお、本実施形態のシミュレーション装置が、上記粉体供給装置等を有する場合、該シミュレーション装置は反応装置等ということもできる。また、粉体供給装置、反応炉、制御装置と、複数の粒子を含む粉体の挙動を解析するための装置とは切り離し可能に構成することもでき、シミュレーション結果を制御装置に反映させた後、切り離して使用しても良い。
［シミュレーション方法］
　次に、本実施形態のシミュレーション方法について説明する。本実施形態のシミュレーション方法は、例えば既述のシミュレーション装置を用いて実施できる。このため、既に説明した事項の一部は説明を省略する。
　本実施形態のシミュレーション方法は、複数の粒子を含む粉体の挙動を解析するためのシミュレーション方法に関する。本実施形態のシミュレーション方法は、図５に示したフロー図に従って実施することができ、以下の工程を有することができる。
　粉体に関連するパラメータを含む第１パラメータを取得する第１パラメータ取得工程（Ｓ１）。　
　複数個の粒子から構成される粒子群を粗視化し、１個の粗視化粒子とした場合の、粗視化粒子についてのパラメータである第２パラメータを算出する第２パラメータ算出工程（Ｓ２）。　
　第１パラメータ、および第２パラメータに基づいて、粗視化粒子の挙動を解析する粗視化粒子挙動解析工程（Ｓ３）。　
　そして、第２パラメータ算出工程（Ｓ２）は、粒子群の弾性エネルギーと、粗視化粒子の弾性エネルギーとの関係を用いた特性方程式の解を用いて、第２パラメータを算出できる。
　各工程について以下に説明する。
（１）第１パラメータ取得工程（Ｓ１）
　第１パラメータ取得工程（Ｓ１）では、解析の対象となる粉体に関連するパラメータを含む第１パラメータを取得できる。既述のシミュレーション装置を用いる場合、例えば第１パラメータ取得部４２１において、第１パラメータ取得工程を実施できる。
　第１パラメータは解析の内容に応じて選択できるため、その具体的な種類は特に限定されない。第１パラメータとしては、離散要素法計算に必要となる各種パラメータが挙げられる。第１パラメータの具体的な例はシミュレーション装置で既に説明したため、ここでは説明を省略する。
　第１パラメータは、データベース等に収録されているデータであってもよく、予め実験を行い求めた実験値であってもよい。また、第１パラメータは、実験結果からシミュレーション等によりフィッティングして算出した計算値であってもよい。
（２）第２パラメータ算出工程（Ｓ２）
　本実施形態のシミュレーション方法では、計算量を抑制するため、粉体が有する複数個の粒子から構成される粒子群を、１個の粗視化粒子に粗視化し、粒子の数を少なくして計算を行うことができる。
　そこで、第２パラメータ算出工程（Ｓ２）では、「（１）粒子の粗視化、および粗視化粒子の粒子挙動の計算に用いるパラメータについて」の中で説明したように、粗視化前の粒子群の弾性エネルギーと、粗視化粒子の弾性エネルギーとの関係を用いて導出した特性方程式の解であるＫｒを用いて、第２パラメータを算出できる。具体的には例えば、粗視化前の粒子群１１全体の弾性エネルギーと、粗視化粒子全体の弾性エネルギーが等しくなると仮定して、既述の垂直方向のＫｒの特性方程式である式（１３）を導出し、式（１３）から垂直方向のＫｒを算出できる。そして、既述の式（１３）に示した特性方程式の解である垂直方向のＫｒを用いて、第２パラメータを算出できる。また、既述の式（１９）に示した接線方向のＫｒの特性方程式を用いて、接線方向のＫｒを算出でき、係る接線方向のＫｒを用いて第２パラメータを算出することもできる。
　既述のように、Ｋｒは、粗視化粒子の挙動を支配するパラメータであり、Ｋｒを用いることで、粗視化粒子の挙動に関連する各種パラメータを算出できる。
　既述のシミュレーション装置を用いる場合、例えば第２パラメータ算出部４２２において、第２パラメータ算出工程を実施できる。
　後述する粗視化粒子挙動解析工程で用いる第２パラメータの種類は解析の内容に応じて選択できるため、特に限定されない。例えば第２パラメータは、粗視化粒子の熱伝導率を含むこともできる。この場合、第２パラメータ算出工程は、既述の特性方程式の解Ｋｒを用いて、熱伝導率を算出できる。
（３）粗視化粒子挙動解析工程（Ｓ３）
　粗視化粒子挙動解析工程（Ｓ３）では、第１パラメータ取得工程（Ｓ１）で得た第１パラメータ、および第２パラメータ算出工程（Ｓ２）で算出した第２パラメータを用いて、粗視化粒子の挙動を解析できる。具体的には離散要素法を用いて計算し、粗視化粒子の挙動を解析できる。粗視化粒子の挙動を解析することで、粉体の挙動を解析することができる。
　なお、ここでいう挙動とは、粗視化粒子の運動による位置の変化だけではなく、温度変化等の状態変化も含む。
（４）初期設定工程
　本実施形態のシミュレーション方法は、例えばさらに初期設定工程を有することもできる。初期設定工程では、解析対象となる粉体を構成する粒子の位置を初期化するとともに、解析の条件、例えば必要に応じて粉体を配置する領域の温度等を設定できる。なお、例えば粗視化粒子挙動解析工程で粗視化粒子の挙動を解析する際に用いるプログラム等に予め初期条件が設定されている場合や、第１パラメータ取得工程により取得する場合には、初期設定工程は実施しなくてもよい。
（５）出力工程
　本実施形態のシミュレーション方法は、例えばさらに出力工程を有することができる。出力工程では、例えば粗視化粒子挙動解析工程（Ｓ３）で得られたシミュレーション結果を、出力部へ出力できる。出力するシミュレーション結果の内容は特に限定されないが、例えば出力部に粗視化された粒子の位置を時系列で画像として出力し、表示することができる。また、例えば出力部に、粉体の温度分布の時系列変化を画像として出力し、表示することもできる。
　以上に説明した本実施形態のシミュレーション方法によれば、複数個の粒子からなる粒子群を１個の粗視化粒子とすることで、計算量を抑制できる。このため、工場で用いるプラント等のような大きな規模での粉体の挙動についても計算量を抑制し、効率的に計算を行うことができる。
　そして、粗視化粒子のパラメータを既述のパラメータＫｒを用いて計算するため、精度よく計算を行うことができる。
　なお、本実施形態のシミュレーション方法は、シミュレーション結果を用いて、粉体を用いた各種製造工程を実施するために、さらに粉体供給工程、反応工程等を有することもできる。
　粉体供給工程では、シミュレーションにより得られた反応炉内での粉体の挙動の結果に基づいて、粉体供給装置から粉体を反応炉へと供給できる。この際粉体としては、シミュレーション結果に基づいて選択した、所定の物性の粉体を供給できる。
　また、反応工程では、粉体供給工程で反応炉へと供給された粉体を加熱できる。この際、シミュレーション結果に基づいて、所定の加熱条件で、粉体を加熱できる。
　上述のように、シミュレーションにより得られた反応炉内での粉体の挙動の結果に基づいて、粉体供給装置から所望の物性の粉体を反応炉へと供給し、さらに反応炉内で所定の加熱条件で加熱することで、粉体の反応割合を高めることができる。また、加熱条件等を最適化できるため、反応の際に用いるエネルギー量を最適化し、生産性を向上できる。
　なお、本実施形態のシミュレーション方法が、上記粉体供給工程等を実施する場合、該シミュレーション方法は粉体の処理方法等ということもできる。
［プログラム］
　次に、本実施形態のプログラムについて説明する。
　本実施形態のプログラムは、複数の粒子を含む粉体の挙動を解析するためのプログラムに関し、コンピュータを以下の各部として機能させることができる。
　粉体に関連するパラメータを含む第１パラメータを取得する第１パラメータ取得部。
　複数個の粒子から構成される粒子群を粗視化し、１個の粗視化粒子とした場合の、粗視化粒子についてのパラメータである第２パラメータを算出する第２パラメータ算出部。
　第１パラメータ、および第２パラメータに基づいて、粗視化粒子の挙動を解析する粗視化粒子挙動解析部。
　第２パラメータ算出部では、粒子群の弾性エネルギーと、粗視化粒子の弾性エネルギーとの関係を用いた特性方程式の解を用いて、第２パラメータを算出させることができる。
　また、第２パラメータは、粗視化粒子の熱伝導率を含むことができ、この場合、第２パラメータ算出部では、既述の特性方程式の解を用いて、熱伝導率を算出できる。
　本実施形態のプログラムは、例えば既述のシミュレーション装置のＲＡＭやＲＯＭ等の主記憶装置または補助記憶装置の各種記憶媒体に記憶させておくことができる。そして、係るプログラムを読み込ませ、ＣＰＵにより実行することにより、ＲＡＭ等におけるデータの読み出しおよび書き込みを行うと共に、入出力インタフェースおよび表示装置を動作させて実行できる。このため、シミュレーション装置で既に説明した事項については説明を省略する。
　上述した本実施形態のプログラムは、インターネットなどのネットワークに接続されたコンピュータ上に格納し、ネットワーク経由でダウンロードさせることで提供してもよい。また、本実施形態のプログラムをインターネットなどのネットワークを介して提供、配布するように構成してもよい。
　本実施形態のプログラムは、ＣＤ－ＲＯＭ等の光ディスクや、半導体メモリ等の記録媒体に格納した状態で流通等させてもよい。
　以上に説明した本実施形態のプログラムによれば、複数個の粒子からなる粒子群を１個の粗視化粒子とすることで、計算量を抑制できる。このため、工場で用いるプラント等のような大きな規模での粉体の挙動についても計算量を抑制し、効率的に計算を行うことができる。
　そして、粗視化粒子のパラメータを既述のパラメータＫｒを用いて計算するため、精度よく計算を行うことができる。
２．第２の実施形態
［シミュレーション装置］
（１）粒子の粗視化、および粗視化粒子の粒子挙動の計算に用いるパラメータについて
　第２の実施形態では、接線方向の運動方程式に関して、オーバーラップ量を算出する際に、回転については角運動量と回転エネルギーが粗視化前後で一致すると仮定する点がここまで説明した第１の実施形態の場合と異なる。また、回転方向の運動方程式に関して、上述のようにして求めたオーバーラップ量を用いてトルクを算出できる。これにより計算量を抑制したまま、粗視化後の粒子群の挙動についてより精度よく解析可能になる。
（１－１）オーバーラップ量について
　通常、接線方向のオーバーラップ量は接触開始から接触終了までの速度の接線方向の成分ｖｔ（以下の式中ではｖｔの上に矢印で示されたもの）の時間積分を用いて以下のように求められる。ここでｔは時間、ベクトルｒ（以下の式中ではｒの上に矢印で示されたもの）は、粒子群１１を構成する粒子の中心から接触点までのベクトル、ベクトルω（以下の式中ではωの上に矢印で示されたもの）は粒子群１１の回転ベクトルを示す。
　なお、以下の式中のハットｔ（以下の式中ではｔの上にハットを付して示されたもの）は、接線方向オーバーラップの単位ベクトルを示す。
Figure JPOXMLDOC01-appb-M000029
 　このため、粒子と壁面との間での各オーバーラップ量は以下の式（２９）で表される。なお、以下の式中の添え字のｔは接線方向成分を意味する。
Figure JPOXMLDOC01-appb-M000030
 　また、粒子間のオーバーラップ量は以下の式（３０）で表される。
Figure JPOXMLDOC01-appb-M000031
 　一方、粗視化粒子２１のオーバーラップ量は以下の式（３１）で表される。なお、ベクトルωｃw（以下の式中ではωｃｗの上に矢印で示されたもの）は粗視化粒子２１の回転ベクトルである。
Figure JPOXMLDOC01-appb-M000032
 　ここで、垂直方向における運動方程式である式（１１）の場合と同様に、粗視化前の粒子群１１と、粗視化粒子２１とで重心位置が一致しているとすると、以下の式（３２）の関係を満たすことになる。
Figure JPOXMLDOC01-appb-M000033
 　ここで、粒子の回転を考慮する場合、重心移動に対して回転方向の自由度が残る。仮に重心位置が一致している場合であっても、式（３１）と式（３２）が一致しない場合がある。このため、粒子群１１の回転運動と粗視化粒子２１の回転運動を換算する必要がある。そこで、回転に関しては角運動量と回転エネルギーが粗視化前後で一致すると仮定する。この場合、角運動量は以下のような等式が成り立つ。
Figure JPOXMLDOC01-appb-M000034
 さらに回転運動エネルギーは以下のような等式が成り立つ。
Figure JPOXMLDOC01-appb-M000035
 　上記式（３３）、式（３４）中、右辺の第１項、第２項は、それぞれ粗視化前の粒子群の自転運動成分（Ｓｐｉｎ）と、公転運動成分（Ｏｒｂｉｔ）とを意味する。ベクトルωｃｗ、ベクトルωｓ、ベクトルωｏ（それぞれの式（３３）、式（３４）中ではωｃｗ、ωｓ、ωｏの上に矢印で示されたもの）はそれぞれ粗視化粒子の角速度、粗視化前の粒子群の自転運動の角速度、粗視化前の粒子群の公転運動の角速度を意味する。また、Ｉｃｗ、Ｉｓ、Ｉｏは粗視化粒子の慣性モーメント、粗視化前の粒子群の自転運動の慣性モーメント、粗視化前の粒子群の公転運動の慣性モーメントを意味する。
　上記式（３３）、式（３４）を連立して公転運動成分を消去すると、以下の式（３５）になる。
Figure JPOXMLDOC01-appb-M000036
 　これより、粗視化粒子の接線方向のオーバーラップ量を以下のように定義できる。
Figure JPOXMLDOC01-appb-M000037
 　この定義の下で、垂直方向の場合と同様に以下の式（３７）が成り立つ。
Figure JPOXMLDOC01-appb-M000038
 　ここでは計算の簡便化のため式（３３）、式（３４）、式（３５）としたが、粗視化粒子の形状によってはいくつかの形式が考えられる。例えば立方体状とした粗視化粒子の形状因子を考慮して慣性モーメントを算出しようとすると以下の式（３８）～式（４１）となる。
Figure JPOXMLDOC01-appb-M000039
Figure JPOXMLDOC01-appb-M000040
Figure JPOXMLDOC01-appb-M000041
Figure JPOXMLDOC01-appb-M000042
 　このため、角運動量保存とベクトルωｓ、ωｏについて、ωｓ＝ωｏを仮定すると式（３８）より以下の式（４２）、式（４３）となる。
Figure JPOXMLDOC01-appb-M000043
Figure JPOXMLDOC01-appb-M000044
 　また、あるいは同様に回転の運動エネルギーが保存することとベクトルωｓ、ωｏについて、ωｓ＝ωｏを仮定するなどしてベクトルωｓ、ベクトルωｃｗの関係を求めてもよい。計算対象とする現象に合わせてこれらを選択することができる。
　このようにして、粗視化前の粒子群を構成する粒子のオーバーラップ量も求めることができる。
（１－２）回転方向の運動方程式について
　ここまで説明した接線方向のオーバーラップ量から接線方向の力を算出し、トルクを計算できる。各粒子間や、粒子と壁間での転がり摩擦は粒子にかかる垂直抗力と接触半径の積に比例するトルクが発生するとした。それぞれの転がり摩擦は各粒子で発生するので、粗視化粒子全体での転がり摩擦抵抗は以下の式（４４）のように表せる。
　式（４４）中、ベクトルＴｔｏｔ＿ｆｒｉｃが粗視化粒子全体の転がり摩擦抵抗を、ベクトルＴｗ，ｆｒｉｃが粒子群１１の壁面１２あるいは外部粒子との間の転がり摩擦抵抗を、ベクトルＴｐ，ｆｒｉｃが粒子群１１の内部粒子間の転がり摩擦抵抗を表している。各ベクトルは、以下の式中では文字の上に矢印を付して表記している。なお、粒子群１１では粒子－粒子間に作用反作用により逆向きの接触力が働くため、式（４４）に示した様に、ベクトルＴｐ，ｆｒｉｃは２倍となる。
　そして、ベクトルＴｗ，ｆｒｉｃと、ベクトルＴｐ，ｆｒｉｃとは、式（４５）、式（４６）により表される。なお、式中のハットωは、回転方向の単位ベクトルを意味する。
Figure JPOXMLDOC01-appb-M000045
Figure JPOXMLDOC01-appb-M000046
Figure JPOXMLDOC01-appb-M000047
 　式（４５）、式（４６）中の接触半径ｒｗ、ｒｐはオーバーラップ量から幾何学的に算出することができ、既述の接線方向の運動方程式の場合の定義に基づいて算出したオーバーラップ量を用いることで精度を高めることができる。また、μｗは粒子群１１の壁面１２あるいは外部粒子との転がり摩擦係数を、μｐは粒子群１１の内部粒子の転がり摩擦係数をそれぞれ意味する。
　そして、回転方向の運動方程式は以下のよう表記できる。
Figure JPOXMLDOC01-appb-M000048
 （２）シミュレーション装置
　本実施形態のシミュレーション装置においても、第２パラメータ算出部は、粒子群の弾性エネルギーと、粗視化粒子の弾性エネルギーとの関係を用いた特性方程式の解を用いて、第２パラメータを算出できる。ただし、接線方向の運動方程式を用いる場合において、オーバーラップ量に関して、回転については角運動量と回転エネルギーが粗視化前後で一致すると仮定し算出できる。また、回転方向の運動方程式を用いる場合において、トルクの算出を行う際に、上記接線方向の運動方程式を用いる場合に算出したオーバーラップ量を用いることができる。
　以上の点以外は第１の実施形態のシミュレーション装置の場合と同様に構成できるため、ここでは説明を省略する。
［シミュレーション方法］
　本実施形態のシミュレーション方法においても、第２パラメータ算出工程では、粗視化前の粒子群の弾性エネルギーと、粗視化粒子の弾性エネルギーとの関係を用いて導出した特性方程式の解を用いて、第２パラメータを算出できる。ただし、接線方向の運動方程式を用いる場合において、オーバーラップ量に関して、回転については角運動量と回転エネルギーが粗視化前後で一致すると仮定し算出できる。また、回転方向の運動方程式を用いる場合において、トルクの算出を行う際に、上記接線方向の運動方程式を用いる場合に算出したオーバーラップ量を用いることができる。
　以上の点以外は第１の実施形態のシミュレーション方法の場合と同様に構成できるため、ここでは説明を省略する。
［プログラム］
　本実施形態のプログラムにおいても、第２パラメータ算出部では、粗視化前の粒子群の弾性エネルギーと、粗視化粒子の弾性エネルギーとの関係を用いて導出した特性方程式の解を用いて、第２パラメータを算出できる。ただし、接線方向の運動方程式を用いる場合において、オーバーラップ量に関して、回転については角運動量と回転エネルギーが粗視化前後で一致すると仮定し算出できる。また、回転方向の運動方程式を用いる場合において、トルクの算出を行う際に、上記接線方向の運動方程式を用いる場合に算出したオーバーラップ量を用いることができる。
　以上の点以外は第１の実施形態のプログラムの場合と同様に構成できるため、ここでは説明を省略する。
　以下に具体的な実施例を挙げて説明するが、本発明はこれらの実施例に限定されるものではない。
［実験例１］
［比較例１－１］
　直方体の容器内に充填した粉体層の温度変化について、表１に示したパラメータを用いて、離散要素法計算により底部から加熱を行った場合の解析を行った。解析により求めた容器内の粉体層の平均温度の変化を図６に示す。なお、反発係数を０．１、摩擦係数を０．７、転がり摩擦係数を０．００１とした。
Figure JPOXMLDOC01-appb-T000049
 ［比較例１－２］
　粗視化倍率αを２とし、８個の粒子から構成される粒子群を１個の粗視化粒子とした点以外は比較例１－１と同様の条件で、直方体の容器内に充填した粉体層の温度変化について解析を行った。なお、粗視化粒子としたため、表１に示すように、粒子の粒径を変更しているが、粒径以外は比較例１－１と同じパラメータを用いて、解析を行った。解析により求めた容器内の粉体層の平均温度の変化を図６に示す。
［実施例１－１］
　粗視化倍率αを２とし、既述のシミュレーション装置を用いて、直方体の容器内に充填した粉体層の温度変化について解析を行った。
　具体的には、第１パラメータ取得部４２１により、解析の対象となる粉体に関連するパラメータを含む第１パラメータを取得した（第１パラメータ取得工程：Ｓ１）。
　この際には、比較例１－２と同じパラメータを取得した。
　次いで、第２パラメータ算出部４２２により、粗視化粒子についてのパラメータである第２パラメータを算出した（第２パラメータ算出工程：Ｓ２）。
　この際、既述の式（１３）、式（１９）に示した特性方程式の解であるパラメータＫｒを用いた式（２６）、式（２７）により、粗視化粒子に関する熱伝導率を算出した。また、特性方程式を用いて反発係数、摩擦係数、転がり摩擦係数を算出した。既述のように、計算に適用するモデルに応じて、これらの係数は調整可能であり、上述のように特性方程式を用いて算出、換算できる。算出したパラメータＫｒを表２に示す。
Figure JPOXMLDOC01-appb-T000050
 　以上の取得、算出した第１パラメータ、および第２パラメータを表１に示す。そして、表１に示したパラメータに基づいて、粗視化粒子の挙動、具体的には温度変化を粗視化粒子挙動解析部４２３において解析した（粗視化粒子挙動解析工程：Ｓ３）。結果を図６に示す。
　図６に示した結果によると、実施例１－１の結果は、粗視化を行っていない、比較例１－１の結果とほぼ一致しており、粗視化粒子のパラメータを適切に選択、設定し、解析を行えていることを確認できた。
　実施例１－１では粗視化を行っているため、粒子の数が比較例１－１と比較して１／８倍になっており、比較例１－１と比較して計算量を抑制できている。
［実験例２］
［比較例２－１］
　回転体であるキルン内の粉体の挙動について、表３に示したパラメータを用いて、離散要素法計算により解析を行った。なお、反発係数を０．７５、摩擦係数を０．３、転がり摩擦係数を０．５とした。
解析により求めたキルン内の粉体の動きを図８Ａ～図８Ｄに、解析により求めたキルン内の粉体の温度分布を図１１Ａ～図１１Ｃに、解析により求めたキルン内の粉体（粒子）の平均温度を図１３にそれぞれ示す。
　図８Ａはキルン７０の回転を開始する前の状態を示しており、グループ分けした第１粉体群７１と、第２粉体群７２とが半分ずつ入っている。図８Ｂ、図８Ｃ、図８Ｄは、それぞれキルンの回転開始後２秒経過時、４秒経過時、６秒経過時の状態を示しており、第１粉体群７１と、第２粉体群７２とが混合されている状態を示している。
　図１１Ａはキルンの回転を開始する前の状態を示しており、温度が第１温度域１０１内にあり、均一であることを確認できる。図１１Ｂ、図１１Ｃは、それぞれキルンの回転開始後３秒経過時、６秒経過時の状態を示しており、キルン７０の外壁側から加熱を行っているため、粉体の中心部側から順に第１温度域１０１、第２温度域１０２、第３温度域１０３が分布していることを確認できる。なお、第１温度域１０１、第２温度域１０２、第３温度域１０３の順に温度が高くなる。
Figure JPOXMLDOC01-appb-T000051
 ［比較例２－２］
　粗視化倍率αを４とし、６４個の粒子から構成される粒子群を１個の粗視化粒子とした点以外は比較例２－１と同様の条件で、回転体であるキルン内の粉体の挙動について解析を行った。なお、表３に示すように、粒径以外は比較例２－１と同じパラメータを用いて、解析を行った。解析により求めたキルン内の粉体の動きを図９Ａ～図９Ｄに、解析により求めたキルン内の粉体の温度分布を図１２Ａ～図１２Ｃに、解析により求めたキルン内の粉体の平均温度を図１３にそれぞれ示す。
　図９Ａはキルン７０の回転を開始する前の状態を示しており、グループ分けした第１粉体群７１と、第２粉体群７２とが半分ずつ入っている。図９Ｂ、図９Ｃ、図９Ｄは、それぞれキルンの回転開始後２秒経過時、４秒経過時、６秒経過時の状態を示しており、第１粉体群７１と、第２粉体群７２とが混合されている状態を示している。
　図１２Ａはキルンの回転を開始する前の状態を示しており、温度が第１温度域１０１内にあり、均一であることを確認できる。図１２Ｂ、図１２Ｃは、それぞれキルンの回転開始後３秒経過時、６秒経過時の状態を示しており、キルン７０の外壁側から加熱を行っているため、粉体の中心部側から順に第１温度域１０１、第２温度域１０２、第３温度域１０３が分布していることを確認できる。なお、第１温度域１０１、第２温度域１０２、第３温度域１０３の順に温度が高くなる。
［実施例２－１］
　粗視化倍率αを４とし、既述のシミュレーション装置を用いて、回転体であるキルン内の粉体の挙動について解析を行った。
　具体的には、第１パラメータ取得部４２１により、解析の対象となる粉体に関連するパラメータを含む第１パラメータを取得した（第１パラメータ取得工程：Ｓ１）。
　この際には、比較例２－２と同じパラメータを取得した。
　次いで、第２パラメータ算出部４２２により、粗視化粒子についてのパラメータである第２パラメータを算出した（第２パラメータ算出工程：Ｓ２）。
　この際、既述の式（１３）、式（１９）に示した特性方程式の解であるパラメータＫｒを用いた式（２６）、式（２７）により、粗視化粒子に関する熱伝導率を算出した。また、特性方程式を用いて反発係数、摩擦係数、転がり摩擦係数を算出した。算出したパラメータＫｒを表４に示す。
Figure JPOXMLDOC01-appb-T000052
　以上のようにして取得、算出した、第１パラメータ、および第２パラメータを表３に示す。そして、表３に示したパラメータに基づいて、粗視化粒子の挙動、具体的にはキルン内の動きや、温度変化を粗視化粒子挙動解析部４２３において解析した（粗視化粒子挙動解析工程：Ｓ３）。解析により求めたキルン内の粉体の動きを図７Ａ～図７Ｄに、解析により求めたキルン内の粉体の温度分布を図１０Ａ～図１０Ｃに、解析により求めたキルン内の粉体の平均温度を図１３にそれぞれ示す。
　図７Ａはキルン７０の回転を開始する前の状態を示しており、グループ分けした第１粉体群７１と、第２粉体群７２とが半分ずつ入っている。図７Ｂ、図７Ｃ、図７Ｄは、それぞれキルンの回転開始後２秒経過時、４秒経過時、６秒経過時の状態を示しており、第１粉体群７１と、第２粉体群７２とが混合されている状態を示している。
　図１０Ａはキルンの回転を開始する前の状態を示しており、温度が第１温度域１０１内にあり、均一であることを確認できる。図１０Ｂ、図１０Ｃは、それぞれキルンの回転開始後３秒経過時、６秒経過時の状態を示しており、キルン７０の外壁側から加熱を行っているため、粉体の中心部側から順に第１温度域１０１、第２温度域１０２、第３温度域１０３が分布していることを確認できる。なお、第１温度域１０１、第２温度域１０２、第３温度域１０３の順に温度が高くなる。
　図７Ａ～図７Ｄと図８Ａ～図８Ｄの比較や、図１０Ａ～図１０Ｃと図１１Ａ～図１１Ｃの比較から明らかなように、実施例２－１の結果は、粗視化を行っていない、比較例２－１の結果とほぼ一致しており、粗視化粒子のパラメータを適切に選択、設定し、解析を行えていることを確認できた。
　図１３においても実施例２－１と比較例２－１とは変曲点を有するなど同様の傾向を示すことが確認できた。
　そして、実施例２－１では粗視化を行っているため、粒子の数が比較例２－１と比較して１／６４倍になっており、比較例２－１と比較して計算量を抑制できている。
　以上に、シミュレーション装置、シミュレーション方法、およびプログラムを、実施形態および実施例等で説明したが、本発明は上記実施形態および実施例等に限定されない。特許請求の範囲に記載された本発明の要旨の範囲内において、種々の変形、変更が可能である。
　本出願は、２０２０年６月１日に日本国特許庁に出願された特願２０２０－９５６２０号、および２０２０年９月２８日に日本国特許庁に出願された特願２０２０－１６２７１９号に基づく優先権を主張するものであり、特願２０２０－９５６２０号、および特願２０２０－１６２７１９号の全内容を本国際出願に援用する。
１１　　　　　　　　　　粒子群
１１Ａ、１１Ｂ　　　　　粒子
２１　　　　　　　　　　粗視化粒子
３０　　　　　　　　　　シミュレーション装置
４２１　　　　　　　　　第１パラメータ取得部
４２２　　　　　　　　　第２パラメータ算出部
４２３　　　　　　　　　粗視化粒子挙動解析部
Ｓ１　　　　　　　　　　第１パラメータ取得工程
Ｓ２　　　　　　　　　　第２パラメータ算出工程
Ｓ３　　　　　　　　　　粗視化粒子挙動解析工程

Claims (7)

  1. 　複数の粒子を含む粉体の挙動を解析するためのシミュレーション装置であって、
    　前記粉体に関連するパラメータを含む第１パラメータを取得する第１パラメータ取得部と、
    　複数個の前記粒子から構成される粒子群を粗視化し、１個の粗視化粒子とした場合の、前記粗視化粒子についてのパラメータである第２パラメータを算出する第２パラメータ算出部と、
    　前記第１パラメータ、および前記第２パラメータに基づいて、前記粗視化粒子の挙動を解析する粗視化粒子挙動解析部と、を有し、
    　前記第２パラメータ算出部は、前記粒子群の弾性エネルギーと、前記粗視化粒子の弾性エネルギーとの関係を用いた特性方程式の解を用いて、前記第２パラメータを算出するシミュレーション装置。
  2. 　前記第２パラメータは、前記粗視化粒子の熱伝導率を含んでおり、
    　前記第２パラメータ算出部は、前記特性方程式の解を用いて、前記熱伝導率を算出する請求項１に記載のシミュレーション装置。
  3. 　回転体内での、前記粉体の挙動を解析する請求項１または請求項２に記載のシミュレーション装置。
  4. 　複数の粒子を含む粉体の挙動を解析するためのシミュレーション方法であって、
    　前記粉体に関連するパラメータを含む第１パラメータを取得する第１パラメータ取得工程と、
    　複数個の前記粒子から構成される粒子群を粗視化し、１個の粗視化粒子とした場合の、前記粗視化粒子についてのパラメータである第２パラメータを算出する第２パラメータ算出工程と、
    　前記第１パラメータ、および前記第２パラメータに基づいて、前記粗視化粒子の挙動を解析する粗視化粒子挙動解析工程と、を有し、
    　前記第２パラメータ算出工程は、前記粒子群の弾性エネルギーと、前記粗視化粒子の弾性エネルギーとの関係を用いた特性方程式の解を用いて、前記第２パラメータを算出するシミュレーション方法。
  5. 　前記第２パラメータは、前記粗視化粒子の熱伝導率を含んでおり、
    　前記第２パラメータ算出工程では、前記特性方程式の解を用いて、前記熱伝導率を算出する請求項４に記載のシミュレーション方法。
  6. 　複数の粒子を含む粉体の挙動を解析するためのプログラムであって、
    　コンピュータを、
    　前記粉体に関連するパラメータを含む第１パラメータを取得する第１パラメータ取得部と、
    　複数個の前記粒子から構成される粒子群を粗視化し、１個の粗視化粒子とした場合の、前記粗視化粒子についてのパラメータである第２パラメータを算出する第２パラメータ算出部と、
    　前記第１パラメータ、および前記第２パラメータに基づいて、前記粗視化粒子の挙動を解析する粗視化粒子挙動解析部と、して機能させ、
    　前記第２パラメータ算出部では、前記粒子群の弾性エネルギーと、前記粗視化粒子の弾性エネルギーとの関係を用いた特性方程式の解を用いて、前記第２パラメータを算出させるプログラム。
  7. 　前記第２パラメータは、前記粗視化粒子の熱伝導率を含んでおり、
    　前記第２パラメータ算出部は、前記特性方程式の解を用いて、前記熱伝導率を算出する請求項６に記載のプログラム。
PCT/JP2021/020745 2020-06-01 2021-05-31 シミュレーション装置、シミュレーション方法、プログラム Ceased WO2021246378A1 (ja)

Priority Applications (4)

Application Number Priority Date Filing Date Title
EP21818068.5A EP4160184B1 (en) 2020-06-01 2021-05-31 Simulation device, simulation method, and program
CN202180039369.1A CN115698672B (zh) 2020-06-01 2021-05-31 模拟装置、模拟方法及存储介质
KR1020227041975A KR102660215B1 (ko) 2020-06-01 2021-05-31 시뮬레이션 장치, 시뮬레이션 방법, 프로그램
US18/058,808 US11892388B2 (en) 2020-06-01 2022-11-25 Simulation device, simulation method, and program

Applications Claiming Priority (4)

Application Number Priority Date Filing Date Title
JP2020-095620 2020-06-01
JP2020095620 2020-06-01
JP2020162719A JP7101387B2 (ja) 2020-06-01 2020-09-28 シミュレーション装置、シミュレーション方法、プログラム
JP2020-162719 2020-09-28

Related Child Applications (1)

Application Number Title Priority Date Filing Date
US18/058,808 Continuation-In-Part US11892388B2 (en) 2020-06-01 2022-11-25 Simulation device, simulation method, and program

Publications (1)

Publication Number Publication Date
WO2021246378A1 true WO2021246378A1 (ja) 2021-12-09

Family

ID=78830357

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/JP2021/020745 Ceased WO2021246378A1 (ja) 2020-06-01 2021-05-31 シミュレーション装置、シミュレーション方法、プログラム

Country Status (5)

Country Link
US (1) US11892388B2 (ja)
EP (1) EP4160184B1 (ja)
KR (1) KR102660215B1 (ja)
CN (1) CN115698672B (ja)
WO (1) WO2021246378A1 (ja)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US11892388B2 (en) * 2020-06-01 2024-02-06 Sumitomo Metal Mining Co., Ltd. Simulation device, simulation method, and program

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN119129075B (zh) * 2024-09-20 2025-12-19 中国矿业大学(北京) 一种用于分析锚杆拉伸承载特性的计算方法

Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2000322407A (ja) * 1999-05-12 2000-11-24 Toyota Central Res & Dev Lab Inc 流体中における媒体の運動解析方法
JP2013210918A (ja) * 2012-03-30 2013-10-10 Sumitomo Chemical Co Ltd 粒子挙動シミュレーション装置、粒子挙動シミュレーション方法、制御プログラム、および記録媒体
CN104636538A (zh) * 2014-12-30 2015-05-20 广东工业大学 基于dem的球磨机颗粒轨迹分析与能耗建模方法
JP2020057135A (ja) 2018-10-01 2020-04-09 住友重機械工業株式会社 シミュレーション装置、シミュレーション方法、及びプログラム
JP2020064450A (ja) * 2018-10-17 2020-04-23 株式会社豊田中央研究所 繊維挙動計算装置、方法、及びプログラム
JP2020095620A (ja) 2018-12-14 2020-06-18 株式会社デンソーテン 画像処理装置および画像処理方法
JP2020162719A (ja) 2019-03-28 2020-10-08 株式会社サンセイアールアンドディ 遊技機
JP2020181393A (ja) * 2019-04-25 2020-11-05 株式会社日立製作所 粒子挙動シミュレーション方法、及び粒子挙動シミュレーションシステム

Family Cites Families (29)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2005529604A (ja) * 2002-06-13 2005-10-06 ザ・トラスティーズ・オブ・ザ・ユニバーシティ・オブ・ペンシルベニア 粗粒度モデルを用いる生体膜をシミュレートするための方法、システム、およびコンピュータプログラム製品
US7813903B2 (en) * 2005-04-13 2010-10-12 Autodesk, Inc. Fixed time step dynamical solver for interacting particle systems
CA2580998A1 (en) * 2006-03-03 2007-09-03 Queen's University At Kingston Adaptive analysis methods
US7565276B2 (en) * 2006-04-05 2009-07-21 Seoul National University Industry Foundation Method of simulating detailed movements of fluids using derivative particles
KR100984048B1 (ko) * 2008-10-14 2010-09-30 한국전자통신연구원 파티클 유체 시뮬레이션에서의 강성체 상호작용 처리 방법
JP5527572B2 (ja) * 2008-10-29 2014-06-18 株式会社Ｉｈｉ シミュレーション方法及びシミュレーション装置
JP2010108183A (ja) * 2008-10-29 2010-05-13 Ihi Corp シミュレーション方法及びシミュレーション装置
JP5241468B2 (ja) * 2008-12-19 2013-07-17 住友重機械工業株式会社 シミュレーション方法及びプログラム
US20100185420A1 (en) * 2009-01-18 2010-07-22 Ejiang Ding Computer system for computing the motion of solid particles in fluid
EP2484648A4 (en) * 2009-09-28 2014-01-15 Asahi Glass Co Ltd LAMINATED GLASS SUBSTRATE, METHOD FOR PRODUCING THE LAMINATED GLASS SUBSTRATE AND ELECTRONIC DEVICE EQUIPPED WITH THE LAMINATED GLASS SUBSTRATE
JP2011081530A (ja) * 2009-10-06 2011-04-21 Univ Of Tokyo 粒子モデルの離散要素法解析シミュレーション方法、離散要素法解析シミュレーションプログラム、及び離散要素法解析シミュレーション装置
US20120150496A1 (en) * 2010-12-09 2012-06-14 Jiun-Der Yu Simplified Fast 3D Finite Difference Particulate Flow Algorithms for Liquid Toner Mobility Simulations
JP5669589B2 (ja) * 2011-01-19 2015-02-12 住友重機械工業株式会社 解析装置
US9183676B2 (en) * 2012-04-27 2015-11-10 Microsoft Technology Licensing, Llc Displaying a collision between real and virtual objects
JP6186675B2 (ja) * 2012-07-05 2017-08-30 日本電気硝子株式会社 ガラス樹脂積層体
JP5930987B2 (ja) * 2013-02-27 2016-06-08 住友重機械工業株式会社 解析装置およびコンピュータプログラム
US9824192B2 (en) * 2013-04-04 2017-11-21 Sumitomo Rubber Industries, Ltd. Simulation method for macromolecular material
US20140358505A1 (en) * 2013-05-31 2014-12-04 The Board Of Trustees Of The University Of Illinois Collision impulse derived discrete element contact force determination engine, method, software and system
US10210287B2 (en) * 2013-12-31 2019-02-19 Disney Enterprises, Inc. Augmented material point method for simulating phase changes and varied materials
US10019827B2 (en) * 2013-12-31 2018-07-10 Disney Enterprises, Inc. Material point method for simulation of granular materials
US10319132B2 (en) * 2014-05-08 2019-06-11 Nvidia Corporation Method and system for representing objects with velocity-dependent particles
US10423736B2 (en) * 2015-08-28 2019-09-24 University Of British Columbia Methods and systems for simulating hydrodynamics in gas-solid fluidized beds
CN107341277A (zh) * 2016-04-29 2017-11-10 中国科学院微电子研究所 基于粒子群优化算法的参数拟合方法及参数拟合系统
JP2018018164A (ja) * 2016-07-25 2018-02-01 富士通株式会社 粒子シミュレーションプログラム、粒子シミュレーション方法、及び情報処理装置
US11062060B2 (en) * 2016-09-09 2021-07-13 Disney Enterprises, Inc. Guided simulation with generalized collision objects
CN107423498B (zh) * 2017-07-13 2020-03-10 山东大学 一种高致密度离散颗粒多相体系的建模方法
JPWO2019181541A1 (ja) * 2018-03-19 2021-03-18 住友重機械工業株式会社 シミュレーション方法、シミュレーション装置、及びプログラム
CN110750833A (zh) * 2019-03-22 2020-02-04 大连理工大学 一种基于改进的移动粒子半隐式法和模态叠加方法求解强非线性时域水弹性问题设计方法
CN115698672B (zh) * 2020-06-01 2024-05-31 住友金属矿山株式会社 模拟装置、模拟方法及存储介质

Patent Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2000322407A (ja) * 1999-05-12 2000-11-24 Toyota Central Res & Dev Lab Inc 流体中における媒体の運動解析方法
JP2013210918A (ja) * 2012-03-30 2013-10-10 Sumitomo Chemical Co Ltd 粒子挙動シミュレーション装置、粒子挙動シミュレーション方法、制御プログラム、および記録媒体
CN104636538A (zh) * 2014-12-30 2015-05-20 广东工业大学 基于dem的球磨机颗粒轨迹分析与能耗建模方法
JP2020057135A (ja) 2018-10-01 2020-04-09 住友重機械工業株式会社 シミュレーション装置、シミュレーション方法、及びプログラム
JP2020064450A (ja) * 2018-10-17 2020-04-23 株式会社豊田中央研究所 繊維挙動計算装置、方法、及びプログラム
JP2020095620A (ja) 2018-12-14 2020-06-18 株式会社デンソーテン 画像処理装置および画像処理方法
JP2020162719A (ja) 2019-03-28 2020-10-08 株式会社サンセイアールアンドディ 遊技機
JP2020181393A (ja) * 2019-04-25 2020-11-05 株式会社日立製作所 粒子挙動シミュレーション方法、及び粒子挙動シミュレーションシステム

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
PATIL, A. V.: "Comparison of CFD-DEM heat transfer simulations with infrared/visual measurements", CHEMICAL ENGINEERING JOURNAL, vol. 277, 2015, pages 388 - 401, XP029611011, DOI: 10.1016/j.cej.2015.04.131 *
See also references of EP4160184A4
YAMANOI MIKIO, NAKATA YOICHI: "DEM Simulation Based on Scaling Rule", JOURNAL OF THE SOCIETY OF POWDER TECHNOLOGY, vol. 55, no. 2, 1 January 2018 (2018-01-01), Kyoto , pages 95 - 103, XP055881921, ISSN: 0386-6157, DOI: 10.4164/sptj.55.95 *

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US11892388B2 (en) * 2020-06-01 2024-02-06 Sumitomo Metal Mining Co., Ltd. Simulation device, simulation method, and program

Also Published As

Publication number Publication date
US20230091287A1 (en) 2023-03-23
CN115698672A (zh) 2023-02-03
EP4160184A1 (en) 2023-04-05
KR102660215B1 (ko) 2024-04-23
EP4160184A4 (en) 2024-06-26
KR20220163531A (ko) 2022-12-09
EP4160184B1 (en) 2025-07-02
US11892388B2 (en) 2024-02-06
CN115698672B (zh) 2024-05-31

Similar Documents

Publication Publication Date Title
Ozel et al. Towards filtered drag force model for non-cohesive and cohesive particle-gas flows
Saeed et al. Mixing study of non-spherical particles using DEM
Elskamp et al. A strategy to determine DEM parameters for spherical and non-spherical particles
Delong et al. Brownian dynamics of confined rigid bodies
(Yali) Tang et al. A new drag correlation from fully resolved simulations of flow past monodisperse static arrays of spheres
Liu et al. CFD-DEM simulation of liquid-solid fluidized bed with dynamic restitution coefficient
Kemp et al. Modelling of particulate drying in theory and practice
US11593541B2 (en) Particle behavior simulation method and particle behavior simulation system
WO2021246378A1 (ja) シミュレーション装置、シミュレーション方法、プログラム
JP7101387B2 (ja) シミュレーション装置、シミュレーション方法、プログラム
Pezo et al. DEM/CFD approach for modeling granular flow in the revolving static mixer
Xie et al. Numerical realization and structure characterization on random close packings of cuboid particles with different aspect ratios
Zhu et al. Numerical investigation of steady and unsteady state hopper flows
Siegmann et al. Efficient Discrete Element Method Simulation Strategy for Analyzing Large‐Scale Agitated Powder Mixers
Sinha et al. A novel computational approach coupled with machine learning to predict the extent of agglomeration in particulate processes
Taghizadeh et al. Numerical analysis of restitution coefficient, rotational speed and particle size effects on the hydrodynamics of particles in a rotating drum: A. Taghizadeh et al.
Shi et al. Numerical study on the elucidation of powder mixing mechanism in a container blender
Girnth et al. Dimensionless quantities in discrete element method: powder model parameterization for additive manufacturing
CN115345088B (zh) 多相流耦合计算方法、装置、设备及存储介质
JP5484127B2 (ja) 粒子挙動解析装置、及び粒子挙動解析方法
JP5527572B2 (ja) シミュレーション方法及びシミュレーション装置
WO2023058764A1 (ja) シミュレーション装置、シミュレーション方法、プログラム
JP2002109445A (ja) 容器内粒子の挙動シミュレーション方法及び挙動シミュレーション装置並びにコンピュータ読み取り可能な記録媒体
JP5669589B2 (ja) 解析装置
Chilamkurti Towards Understanding the Heat Transfer Behavior of Dense Granular Media

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 21818068

Country of ref document: EP

Kind code of ref document: A1

ENP Entry into the national phase

Ref document number: 20227041975

Country of ref document: KR

Kind code of ref document: A

NENP Non-entry into the national phase

Ref country code: DE

ENP Entry into the national phase

Ref document number: 2021818068

Country of ref document: EP

Effective date: 20230102

WWG Wipo information: grant in national office

Ref document number: 2021818068

Country of ref document: EP