WO2017150626A1 - 粒子シミュレーション装置、粒子シミュレーション方法及び粒子シミュレーションプログラム - Google Patents

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

Info

Publication number
WO2017150626A1
WO2017150626A1 PCT/JP2017/008168 JP2017008168W WO2017150626A1 WO 2017150626 A1 WO2017150626 A1 WO 2017150626A1 JP 2017008168 W JP2017008168 W JP 2017008168W WO 2017150626 A1 WO2017150626 A1 WO 2017150626A1
Authority
WO
WIPO (PCT)
Prior art keywords
particle
distance
virtual
interaction force
calculating
Prior art date
Application number
PCT/JP2017/008168
Other languages
English (en)
French (fr)
Inventor
泰介 西浦
Original Assignee
国立研究開発法人海洋研究開発機構
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by 国立研究開発法人海洋研究開発機構 filed Critical 国立研究開発法人海洋研究開発機構
Priority to CN201780014249.XA priority Critical patent/CN108701275A/zh
Priority to AU2017227323A priority patent/AU2017227323B2/en
Priority to GB1815732.1A priority patent/GB2564326B/en
Priority to US16/080,875 priority patent/US20200110911A1/en
Publication of WO2017150626A1 publication Critical patent/WO2017150626A1/ja

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/25—Design optimisation, verification or simulation using particle-based methods
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/10—Numerical modelling

Definitions

  • the present invention relates to a particle simulation apparatus, a particle simulation method, and a particle simulation program for simulating the behavior of particles and structures.
  • the distance function is handled as one of the following, for example, with respect to the movement of the structure with time.
  • the distance function is fixed, and the coordinate system of the particle is translated and rotated corresponding to the movement and rotation of the structure.
  • this method is not applicable when there are a plurality of structures, that is, when there are a plurality of moving axes and rotating axes.
  • (2) Recreate the distance function each time the structure moves.
  • this method requires enormous calculation cost for recreating the distance function, and is not easy to use for actual calculation.
  • the advection calculation of the distance function using a differential equation is performed. However, since this method is a lattice method, numerical diffusion occurs and the distance function is lost.
  • the present invention has been made in view of the above problems, and accurately calculates the distance between the particle and the structure with a small calculation load, thereby appropriately performing the simulation of the particle and the structure.
  • An object is to provide a particle simulation apparatus, a particle simulation method, and a particle simulation program.
  • a particle simulation apparatus calculates the position and velocity of a particle based on the interaction force between the particle and the structure in the work space.
  • a particle simulation apparatus for simulating the behavior of the particles and the structure the structure initial information acquiring means for acquiring the structure initial information indicating the initial position and shape of the structure, and the structure initial information acquiring means Based on the initial position and shape of the structure indicated by the acquired structure initial information, a plurality of virtual areas are set in the vicinity of the surface of the structure, and for each of the set virtual areas, up to the structure Virtual area setting means for calculating the distance, position information acquisition means for acquiring position information indicating the position of the particles, and position information acquired by the position information acquisition means Based on the distance calculated by the virtual area setting means for the virtual area specifying means for specifying the virtual area within a preset distance from the position of the particle indicated by the virtual area specifying means, A particle distance calculating means for calculating a distance from the particle to the structure, and an interaction force calculating
  • the distance between the particle and the structure is calculated based on the distance calculated for a plurality of virtual regions set in the vicinity of the surface of the structure.
  • the virtual area is moved while maintaining the positional relationship with the structure. Therefore, it is not necessary to recalculate the distance (distance function) for the virtual region.
  • the distance only needs to be calculated only for the virtual region set in the vicinity of the structure.
  • the virtual region is set for each structure, it can be applied even when there are a plurality of structures, that is, when there are a plurality of different movement axes and rotation axes.
  • the distance between the particle and the structure is accurately calculated with a small calculation load, thereby appropriately performing the simulation of the particle and the structure. be able to.
  • the virtual area setting means may set a plurality of virtual areas in the vicinity of the surface of the structure within a range corresponding to a preset distance. According to this configuration, it is possible to set as many virtual regions as necessary, and it is possible to more efficiently simulate particles and structures.
  • the time transition means may calculate the position and speed of the structure at the next time step based on the interaction force calculated by the interaction force calculation means. According to this configuration, since the influence of the particles on the structure can also be simulated, more accurate simulation of the particles and the structure can be performed.
  • the virtual area may be a particulate area. According to this configuration, the virtual region can be handled appropriately and easily, and the simulation of the particle and the structure can be performed appropriately and easily.
  • the present invention can be described as an invention of a particle simulation apparatus as described above, and can also be described as an invention of a particle simulation method and a particle simulation program as follows. This is substantially the same invention only in different categories, and has the same operations and effects.
  • the particle simulation method calculates the position and velocity of the particle based on the interaction force between the particle and the structure in the work space, and the particle and the structure.
  • a particle simulation method that is an operation method of a particle simulation apparatus that simulates the behavior of a structure, the structure initial information acquisition step for acquiring structure initial information indicating the initial position and shape of the structure, and the structure initial information acquisition step Based on the initial position and shape of the structure indicated by the structure initial information acquired in step 2, a plurality of virtual areas are set in the vicinity of the surface of the structure, and each of the plurality of virtual areas set up to the structure
  • a virtual region setting step for calculating the distance of the position a position information acquisition step for acquiring position information indicating the position of the particles,
  • a virtual region specifying step that specifies a virtual region within a preset distance from the particle position indicated by the position information acquired in the acquiring step, and the virtual region specified in the virtual region specifying step Based on the calculated distance, a particle distance
  • the particle simulation program calculates the position and velocity of the particle based on the interaction force between the particle and the structure in the work space.
  • a particle simulation program for functioning as a particle simulation apparatus for simulating the behavior of the structure wherein the computer obtains structure initial information acquisition means for acquiring structure initial information indicating the initial position and shape of the structure, and the structure Based on the initial position and shape of the structure indicated by the initial structure information acquired by the initial information acquisition means, a plurality of virtual areas are set near the surface of the structure, and each of the set virtual areas is set.
  • a virtual region setting means for calculating the distance to the structure, and position information indicating the position of the particles Position information acquisition means, virtual area specification means for specifying a virtual area within a preset distance from the position of the particle indicated by the position information acquired by the position information acquisition means, and virtual specified by the virtual area specification means For the region, based on the distance calculated by the virtual region setting unit, a particle distance calculation unit that calculates the distance from the particle to the structure, and based on the distance calculated by the particle distance calculation unit, the particle and the structure
  • the interaction force calculation means for calculating the interaction force between the particles, the position and velocity of the particles in the next time step based on the interaction force calculated by the interaction force calculation means, and the structure And function as time transition means for moving the virtual area while maintaining the positional relationship.
  • the distance there is no need to recalculate the distance (distance function) for the virtual region.
  • the distance unlike the conventional distance function, the distance only needs to be calculated only for the virtual region set in the vicinity of the structure.
  • the virtual region since the virtual region is set for each structure, it can be applied even when there are a plurality of structures, that is, when there are a plurality of different movement axes and rotation axes. Further, since the virtual region associated with the distance moves, numerical diffusion due to advection calculation does not occur. That is, according to one embodiment of the present invention, the distance between the particle and the structure can be accurately calculated with a small calculation load, and thereby the simulation of the particle and the structure can be appropriately performed.
  • FIG. 1 shows a particle simulation apparatus 10 according to this embodiment.
  • the particle simulation device 10 is a device that simulates (analyzes) the behavior of a plurality of spherical particles (a shape having a constant distance from the center to the surface) and a structure in a work space (calculation region). Specifically, the particle simulation apparatus calculates the force acting on each particle based on the position and velocity of each particle and structure at each time step that is a simulation time.
  • the force acting on each particle and structure includes, for example, a contact force that is an interaction force due to contact (collision) that is an interaction between particles and the structure and between particles.
  • the particle simulation apparatus 10 calculates the position and velocity of each particle and structure in the next time step based on the calculated force.
  • the particle simulation in this embodiment is performed based on DEM, you may be performed based on other particle methods, such as SPH and MPS. Further, depending on the type of particle method, the particle diameter and rotation speed may not be retained as information.
  • the particles to be simulated by the particle simulation apparatus 10 according to the present embodiment include arbitrary particles that have been targets for conventional particle simulation. For example, earth and sand or powder can be targeted. Or it is good also as a target supposing that fluid and solid consist of a plurality of particles.
  • the structure to be simulated by the particle simulation apparatus 10 according to the present embodiment includes an arbitrary object. For example, an apparatus for agitating and mixing powder (a part for agitating and mixing powder (agitating blade)) is an object. There may be a plurality of structures to be simulated.
  • the physical problem can be simulated by the simulation by the particle simulation apparatus 10 according to the present embodiment.
  • simulation of powder agitation and mixing can be performed.
  • it may be used for verification of drilling capacity according to the shape of the drill bit and for the design of construction equipment and heavy machinery such as a bulldozer.
  • the simulation by the particle simulation apparatus 10 according to the present embodiment can be applied in industrial fields such as civil engineering and powder. It can also be used for the simulation of natural phenomena such as dikes, buildings, terrain and tsunami and earth and sand, landslides and avalanches.
  • the simulation by the particle simulation apparatus 10 according to the present embodiment can be applied in the field of disaster prevention.
  • the particle simulation apparatus 10 is configured as a computer including hardware such as a CPU (Central Processing Unit), a GPU (Graphics Processing Unit), a memory, a hard disk, and a display, for example. When these components operate according to a program or the like, a function as a particle simulation apparatus 10 described later is exhibited.
  • the particle simulation apparatus 10 operates particularly effectively in an apparatus capable of parallel computation. Note that the particle simulation device 10 does not necessarily include a GPU as an arithmetic device, and may have a configuration including only a CPU.
  • the particle simulation apparatus 10 includes a simulation information holding unit 11, a structure initial information acquisition unit 12, a virtual region setting unit 13, a position information acquisition unit 14, a virtual region specification unit 15, and particles.
  • the distance calculation unit 16, the interaction force calculation unit 17, and the time transition unit 18 are configured.
  • the work area which is an area in which particles move, is, for example, a three-dimensional space, and is divided (segmented) into cubic cells (mesh) in which the size of one side is set in advance.
  • the cell may have a shape other than a cube, for example, a rectangular parallelepiped.
  • the particle simulation apparatus 10 previously divides a work space into cells before performing simulation processing, and grasps in advance how the work space is divided into cells.
  • the size of the one side is set in advance according to, for example, the particle diameter of a plurality of particles to be simulated.
  • Each cell in the work area is assigned a cell number that identifies the cell. The cell numbers are given in order according to the position of the cell in the work space, for example.
  • the simulation information holding unit 11 is a means for holding information used for simulation. Specifically, the simulation information holding unit 11 holds particle information for each of a plurality of particles in the work area.
  • the particle information includes information indicating particle coordinates, particle velocity, and particle radius.
  • the coordinate of the particle is a three-dimensional coordinate indicating the position of the particle in the work space.
  • Particle speed includes translational and rotational speeds.
  • information at the start of the simulation is input in advance by the user of the particle simulation apparatus 10 or the like to the simulation information holding unit 11, and information during the simulation is It is updated by the time transition unit 18 described later.
  • the particle radius is previously input to the simulation information holding unit 11 by the user of the particle simulation apparatus 10 or the like.
  • the simulation information holding unit 11 holds structure information for each structure in the work area.
  • the structure information includes information indicating the coordinates of the structure, the speed of the structure, and the shape of the structure.
  • the coordinates of the structure are three-dimensional coordinates indicating the position of the structure in the work space (for example, the position of a specific part of the structure).
  • the speed of the structure includes a translation speed and a rotation speed.
  • the information (initial information) at the start of the simulation is preliminarily stored in the simulation information holding unit 11 by the user of the particle simulation apparatus 10 or the like.
  • the information that has been input and that is being simulated is updated by the time transition unit 18 described later.
  • the information indicating the shape of the structure is, for example, CAD (computer-aided design) data that represents the shape of the structure with polygons.
  • Information indicating the shape of the structure is input in advance to the simulation information holding unit 11 by the user of the particle simulation apparatus 10 or the like.
  • the structure information may include information indicating how the structure moves (for example, a movement axis or a rotation axis). The information is also input in advance to the simulation information holding unit 11 by the user of the particle simulation apparatus 10 or the like.
  • the simulation information holding unit 11 holds virtual particle information as will be described later.
  • maintenance part 11 may input and hold
  • Such information includes a friction coefficient, an elastic coefficient, a viscous damping coefficient, a restitution coefficient, and the like.
  • the structure initial information acquisition unit 12 is a structure initial information acquisition unit that acquires structure initial information indicating the initial position and shape of the structure. For example, the structure initial information acquisition unit 12 acquires initial information related to the structure held in the simulation information holding unit 11 and information indicating the shape of the structure as the structure initial information. The structure initial information acquisition unit 12 outputs the acquired structure initial information to the virtual area setting unit 13.
  • the virtual region setting unit 13 sets a plurality of virtual regions near the surface of the structure based on the initial position and shape of the structure indicated by the structure initial information acquired by the structure initial information acquisition unit 12
  • the virtual area setting means calculates the distance (closest distance) to the structure for each of the set virtual areas.
  • the virtual region is a region used for calculating the distance from the particle to be simulated to the structure.
  • the virtual region is a spherical region (a shape having a constant distance from the center to the surface) and is called a virtual particle.
  • the lengths of the diameters of the plurality of virtual particles are uniform and set in advance.
  • the virtual particles may be larger than the simulation target particles in consideration of the calculation load and the like.
  • the distance from the particle to be simulated to the structure has been obtained using a signed distance function as shown in Non-Patent Document 1.
  • the signed distance function is a function in which a distance from a structure is stored at a lattice point obtained by dividing a work area into cells (by an orthogonal lattice).
  • FIG. 2A shows an isoline diagram in which positions of equal distances from the structure 100 generated based on the lattice point-like distance function values are connected by lines.
  • a plurality of virtual particles 110 are set in the vicinity of the surface of the structure 100 (the boundary between the structure 100 and other portions).
  • the virtual area setting unit 13 first calculates a signed distance function for the structure 100 whose initial position and shape are indicated by the structure initial information. This calculation may be performed using a conventional method as described in Non-Patent Document 2.
  • the virtual area setting unit 13 detects the position of the vicinity from the structure 100 in the work area based on the calculated signed distance function. As shown in FIG. 2B, the vicinity from the structure includes the outside, surface, and inside regions of the structure.
  • the virtual region setting unit 13 sets a plurality of virtual particles 110 so as to come into contact with each other in the neighboring region.
  • the virtual area setting unit 13 sets a lattice point of a simple cubic lattice or a face-centered cubic lattice in the region, and sets (places) the virtual particle 110 at the lattice point.
  • the setting of the virtual particles only needs to know the distance from the structure 100, so that it is not always necessary to use a signed distance function, and may be performed by an arbitrary method.
  • the part of the surface of the structure 100 where the virtual particle 110 is set is a part where the simulation target particle may come into contact (interact).
  • virtual particles 110 are set in the vicinity of the surface over the entire surface of the structure 100 as shown in FIG. .
  • the virtual particles 110 need not be set in the vicinity of the portion.
  • the virtual area setting unit 13 sets a predetermined distance (predetermined thickness) in the vicinity of the structure in the normal direction of the surface of the structure 100. For example, as shown in FIG. 3, a region having a distance H from the surface of the structure 100 to the outside is formed with the distance d as the center, and from there to the outside and inside (vertical direction in FIG. 3). The region is in the vicinity of the object 100.
  • the distance d is a distance between the surface of the structure 100 and the center of the particle 101 where the structure 100 and the particle 101 to be simulated begin to interact.
  • the size of the particle 101 to be simulated is not uniform, for example, for dH in FIG. It is also possible to use the maximum value of.
  • the distance d is the radius of the particle 101.
  • the distance H is a distance for specifying the virtual particle 110 used for calculating the distance from the particle 101 to the structure 100 (details will be described later).
  • the distance H is set in advance.
  • the virtual region setting unit 13 may set a plurality of virtual particles 110 in the vicinity of the surface of the structure 100 in a range corresponding to the distance H.
  • the virtual particles 110 may not be set on the surface and inside of the structure 100.
  • the region (distance and thickness in the normal direction) where the virtual particles 110 are set need not be set (calculated) as described above, and may be set in advance.
  • the positional relationship between the structure 100 and the set virtual particle 110 is fixed throughout the simulation.
  • the virtual region setting unit 13 is a distance from the virtual particle 110 to the structure 100 (related to the surface on which the virtual particle 110 is arranged), for example, the surface of the structure 100 from the center position of the virtual particle 110 The distance to is calculated. For this calculation, any method such as a conventional interpolation method using a signed distance function can be used. The calculation of the distance from the virtual particle 110 to the structure 100 may be performed together with the arrangement of the virtual particle 110.
  • the set virtual particle 110 is uniquely identified by assigning a unique number to each virtual particle 110.
  • the virtual region setting unit 13 stores and holds the information indicating the set position of the virtual particle 110 and the distance in the simulation information holding unit 11.
  • the virtual region setting unit 13 may calculate a normal vector with respect to the surface of the structure 100 when calculating the distance of the virtual particle 110 and associate it with the distance.
  • the information on the normal vector can be used in the interpolation process described later.
  • the virtual particles 110 may be associated with mechanical properties such as Young's modulus and used for processing described later.
  • the processing by the structure initial information acquisition unit 12 and the virtual region setting unit 13 is performed before the calculation related to each time step of the simulation. Subsequent functional units are configured for calculation related to each time step of the simulation.
  • the position information acquisition unit 14 is position information acquisition means for acquiring position information indicating the position of the particle 101 to be simulated. Specifically, the position information acquisition unit 14 acquires the particle information of the current time step of each particle 101 held in the simulation information holding unit 11 as the position information of the particle 101. The position information acquisition unit 14 inputs the acquired particle information of the particles 101 to the virtual region specifying unit 15.
  • the position information acquisition unit 14 also acquires position information indicating the position of the virtual particle 110.
  • the position information acquisition unit 14 acquires information indicating the current time step position of each particle 101 held in the simulation information holding unit 11.
  • the position information acquisition unit 14 outputs the acquired position information of the virtual particles 110 to the virtual region specifying unit 15.
  • the virtual region specifying unit 15 is a virtual region specifying unit that specifies the virtual particle 110 within a preset distance H from the position of the particle 101 indicated by the position information acquired by the position information acquiring unit 14.
  • the virtual particle 110 specified here is used to calculate the distance from the particle 101 to the structure 100 by interpolation. That is, the distance H is a distance indicating a range (interpolation range) in which interpolation is performed, and a value is set in advance so that the distance from the particle 101 to the structure 100 can be appropriately calculated by interpolation.
  • the distance H is set when the virtual particles 110 are set in a face-centered cubic lattice so as to satisfy H> D ⁇ n when the virtual particles 110 are set in a simple cubic lattice. Should satisfy H> D.
  • D is the diameter of the virtual particle 110 as shown in FIG. 3, and n is the number of dimensions of the work area.
  • the virtual region specifying unit 15 determines the position of the particle 101 (particle i) (the position of the center of the particle 101) indicated by the information input from the position information acquisition unit 14 and each virtual A distance P ij between the position of the particle 110 (particle j) (the center position of the virtual particle 110) is calculated.
  • the virtual area specifying unit 15 determines whether or not the calculated distance P ij is within the distance H. Based on the determination, the virtual region specifying unit 15 specifies a virtual particle 110 whose distance from the particle 101 is within the distance H for each particle 101. In the example shown in FIG. 4, four virtual particles 110a are specified.
  • the virtual region specifying unit 15 outputs information indicating the virtual particles 110 specified for each particle 101 to the particle distance calculating unit 16.
  • the identification of the virtual particles 110 within the distance H may be performed efficiently using a method of contact determination between particles (determination of whether the distance between particles is within a threshold) in the conventional particle simulation. For example, based on the cell in which the particle 101 and the virtual particle 110 are located as described in JP-A-2015-115567 (Patent Document 1) and JP-A-2010-238030 (Patent Document 2), the particle 101 And a pair of virtual particles 110 may be set, and a distance may be determined for the pair.
  • the particle distance calculation unit 16 calculates the distance (nearest distance) from the particle 101 to the structure 100 based on the distance calculated by the virtual region setting unit 13 for the virtual particle 110 specified by the virtual region specifying unit 15. It is a particle distance calculation means for calculating.
  • the particle distance calculation unit 16 acquires information indicating the distance to the structure 100 for the virtual particles 110 specified for each particle 101 and held in the simulation information holding unit 11.
  • the particle distance calculation unit 16 determines the distance from the particle 101 to the structure 100 (from the center position of the particle 101 (particle i) shown in FIG. 4 to the surface of the structure 100) based on the distance indicated by the acquired information.
  • the distance dis i is calculated, for example, by the conventional interpolation shown in Non-Patent Document 1, or by any other method, such as SPH.
  • the kernel approximation used can also be used.
  • the particle 101 Does not calculate the distance to the structure 100, assuming that the distance to the structure 100 is not at least the distance at which interaction occurs. In this case, processing using the distance is not performed.
  • the particle distance calculation unit 16 outputs information indicating the distance to the structure 100 calculated for each particle 101 to the interaction force calculation unit 17.
  • the particle distance calculation unit 16 normally calculates a normal vector from the particle 101 to the surface of the structure 100 when the distance is calculated.
  • the interaction force calculated by the interaction force calculation unit 17 can be expressed as a vector quantity.
  • the direction represented by the vector quantity represents the direction in which the interaction force acts.
  • the normal vector can be obtained by calculating a spatial gradient with respect to the distance (distance function data) to the structure 100 for the virtual particle 110 specified for each particle 101.
  • the calculation of the normal vector from the spatial gradient may be performed by a conventional method shown in Non-Patent Document 1, for example.
  • the information may be calculated by interpolation. It may be performed by any other method. For example, a gradient calculation method by kernel approximation used in SPH or the like can also be used.
  • the particle distance calculation unit 16 outputs information indicating the calculated normal vector together with the information indicating the distance to the interaction force calculation unit 17.
  • the interaction force calculator 17 is an interaction force calculator that calculates an interaction force between the particle 101 and the structure 100 based on the distance calculated by the particle distance calculator 16.
  • the interaction force calculation unit 17 determines whether or not the distance from the particle 101 to the structure 100 is smaller than a preset threshold value, and an interaction occurs between the particle 101 and the structure 100. It is determined whether or not to obtain (interaction determination).
  • the interaction force calculation unit 17 determines that the distance is smaller than the threshold value
  • the interaction force calculation unit 17 determines that an interaction occurs between the particle 101 and the structure 100.
  • the interaction force calculation unit 17 determines that the distance is not smaller than the threshold value
  • the interaction force calculation unit 17 determines that no interaction occurs between the particle 101 and the structure 100.
  • the interaction force calculation unit 17 determines whether or not the distance is smaller than the radius d of the particle 101 that is the threshold, and determines that the distance is smaller than the threshold. It is determined that the particle 101 and the structure 100 are in contact with each other and a contact force is generated between them.
  • the calculation of the interaction force can be performed by the same method as in the conventional particle simulation or any other method.
  • the interaction force calculation unit 17 normally calculates the interaction force as a vector amount based on the normal vector from the particle 101 to the surface of the structure 100 calculated by the particle distance calculation unit 16.
  • the interaction force calculation unit 17 may calculate not only the interaction force between the particle 101 and the structure 100 but also the interaction force between the particles 101 and the interaction force between the structures 100. Calculation of the interaction force between the particles 101 can be performed by a method similar to the conventional particle simulation described in Patent Documents 1 and 2 or any other method. Calculation of the interaction force between the structures 100 can also be performed by a method similar to the conventional simulation or any other method. In calculating the interaction force, parameters necessary for calculating the interaction force, which are held in the simulation information holding unit 11, may be acquired and used for calculating the interaction force.
  • the interaction force calculation unit 17 calculates the interaction force between the individual particles 101 and the structure 100 calculated as described above, the interaction force between the particles 101, and the interaction force between the structures 100. The sum of interaction forces for the particle 101 and each structure 100 is calculated. The interaction force calculation unit 17 outputs information indicating the calculated total interaction force of each particle 101 and each structure 100 to the time transition unit 18.
  • the time transition unit 18 calculates the positions and velocities of the particles 101 and the structure 100 in the next time step based on the interaction force calculated by the interaction force calculation unit 17, and the position with the structure 100. This is time transition means for moving the virtual particles 110 while maintaining the relationship. Specifically, the time transition unit 18 acquires the particle information and the structure information held in the simulation information holding unit 11 for each particle 101 and the structure 100, and the current time step indicated in the particle information and the structure information. The coordinates and speed of the next time step are calculated from the coordinates and speed of the current time and the sum of the interaction forces. This calculation can be performed, for example, by a method similar to the conventional simulation or any other method. The time transition unit 18 updates the particle information and structure information held in the simulation information holding unit 11 for each particle 101 and structure 100 at the calculated position and velocity of the particle 101 and structure 100 in the next time step. To do.
  • the time transition unit 18 calculates the position of the next time step for the structure 100
  • the virtual particle 110 set for the structure 100 is maintained in the positional relationship with the structure 100 (the structure Move (similar to the movement of the object 100).
  • the time transition unit 18 updates the information indicating the virtual particle 110 held in the simulation information holding unit 11 with information indicating the position of the moved virtual particle 110 (the information is updated with the virtual particle 110 in the next time step).
  • the information is stored in the simulation information holding unit 11 as information indicating the position).
  • the position and speed of the structure 100 need not be calculated based on the interaction force. For example, when the influence of the particle 101 on the structure 100 is extremely small and no interaction between the structures 100 occurs, the calculation may not be performed. In this case, for example, information indicating structure 100 in advance how it works is, is stored in the simulation information holding unit 11, the time transition 18, the movement (the next time step structures 100 accordingly (Position and speed) may be calculated. Further, the time transition unit 18 moves the virtual particles 110 according to the movement of the structure 100 (according to a predetermined trajectory of the surface of the structure 100).
  • the particle simulation apparatus 10 it is determined whether or not the simulation end condition is satisfied every time the calculation in one time step is completed. For example, when the calculation for the preset number of times (time step) is completed, it is determined that the termination condition is satisfied. When it is determined that the end condition is satisfied, the particle simulation apparatus 10 ends the simulation. In this case, for example, output of calculation results to a display device or another device is performed. If it is determined that the end condition is not satisfied, the calculation of the next time step is repeated.
  • the above is the configuration of the particle simulation apparatus 10.
  • particle simulation method executed by the particle simulation apparatus 10, which is an operation method of the particle simulation apparatus 10 according to the present embodiment, will be described with reference to the flowchart of FIG. This process is started, for example, when the user of the particle simulation apparatus 10 performs an operation for starting the simulation on the particle simulation apparatus 10.
  • the structure initial information acquisition unit 12 acquires structure initial information indicating the initial position and shape of the structure 100 held in the simulation information holding unit 11 (S01, structure initial information). Acquisition step). An example of the shape of the structure 100 is shown in FIG.
  • the acquired structure initial information is output from the structure initial information acquisition unit 12 to the virtual region setting unit 13.
  • the signed distance function of the structure 100 related to the structure initial information is calculated by the virtual area setting unit 13 (S02, virtual area setting step).
  • the virtual region setting unit 13 uses the signed distance function to place (set) a plurality of virtual particles 110 in the vicinity of the surface of the structure 100 (S03, virtual region setting step).
  • An example of the virtual particle 110 arranged in the vicinity of the surface of the structure 100 is shown in FIG.
  • the distance from the virtual particle 110 to the structure 100 is calculated for the virtual particle 110 by the virtual region setting unit 13 (S04, virtual region setting step).
  • Information indicating the position and distance of the virtual particle 110 is held in the simulation information holding unit 11. The above is the process performed before the calculation related to each time step of the simulation.
  • the position information acquisition unit 14 acquires the position information indicating the positions of the particles 101 and the virtual particles 110 to be simulated, which are held in the simulation information holding unit 11 (S05, position information acquisition step).
  • the acquired position information is output from the position information acquiring unit 14 to the virtual area specifying unit 15.
  • the virtual particle 110 within the distance H is specified from the position of the particle 101 indicated by the position information by the virtual region specifying unit 15 (S06, virtual region specifying step).
  • Information indicating the virtual particle 110 specified for each particle 101 is output from the virtual region specifying unit 15 to the particle distance calculating unit 16.
  • the distance from the particle 101 to the structure 100 is calculated by the particle distance calculation unit 16 based on the distance from the particle 101 to the structure 100 for the virtual particle 110 specified for each particle 101 (S07, particle distance). Calculation step). This calculation is performed by interpolation, for example, as described above. Information indicating the distance to the structure 100 calculated for each particle 101 is output from the particle distance calculation unit 16 to the interaction force calculation unit 17.
  • the interaction force calculation unit 17 calculates the interaction force between the particle 101 and the structure 100 based on the distance to the structure 100 calculated for each particle 101 (S08, interaction). Force calculation step). At this time, the interaction force between the particles 101 and the interaction force between the structures 100 may be calculated. Subsequently, the interaction force calculation unit 17 calculates the sum of interaction forces for each particle 101 and each structure 100 from the calculated interaction force. Information indicating the sum of the interaction forces of the calculated particles 101 and the structures 100 is output from the interaction force calculation unit 17 to the time transition unit 18.
  • the time transition unit 18 calculates the positions and velocities of the particle 101 and the structure 100 in the next time step based on the calculated interaction force, and maintains the positional relationship with the structure 100.
  • the virtual particle 110 is moved as it is (S09, time transition step). Based on the above calculation, the particle information, the structure information, and the information indicating the position of the virtual particle 110 held in the simulation information holding unit 11 are updated to the information of the next time step.
  • the particle simulation apparatus 10 it is determined whether or not a simulation end condition is satisfied (S10). If it is determined that the termination condition is satisfied (YES in S10), the process (simulation) is terminated. If it is determined that the termination condition is not satisfied (NO in S10), the time step is advanced by one, and the above-described processing (S05 to S10) in the next time step is performed.
  • FIG. 6 (c) to FIG. 6 (f) show examples of states at each time step.
  • the time step advances from FIG. 6 (c) to FIG. 6 (f).
  • two rotating structures (stirring blades) for stirring the powder are provided inside the container, and the powder falls into the container from above and is stirred by the structure.
  • Each structure can be given different stirring speeds.
  • the above is the processing executed by the particle simulation apparatus 10 according to the present embodiment.
  • the distance between the particle 101 and the structure 100 is calculated based on the distance calculated for the plurality of virtual particles 110 set near the surface of the structure 100. Is called.
  • the virtual particle 110 is moved while maintaining the positional relationship with the structure 100. Therefore, it is not necessary to recalculate the distance (distance function) for the virtual particle 110.
  • each virtual particle 110 is set for each structure 100 and can be moved independently along with the movement of each structure 100. That is, the present invention can be applied even when there are a plurality of structures 100, that is, when there are a plurality of different movement axes and rotation axes.
  • the calculation load is very simple and less than the advection calculation. That is, according to the present embodiment, the distance between the particle 101 and the structure 100 can be accurately calculated with a small calculation load, and thereby the simulation of the particle 101 and the structure 100 can be appropriately performed.
  • the amount of data to be stored (memory) compared to the conventional distance function
  • the amount used can also be reduced primarily.
  • An example of memory usage is shown below.
  • An example in which the size of the work space is 160D ⁇ 200D ⁇ 100D, the interpolation range H is 3.5d, the diameter (lattice size) of the virtual particle 110 is D, and the interaction distance d is 0.5D is shown.
  • the number of particles required for the distance calculation is 85,336, and the total memory usage is 2,730,752 (Byte) (85,336 ⁇ 4 ⁇ 8: a variable of 8 bytes per particle. 4 (x, y, z coordinates and distance function values are required).
  • the number of grids required for the distance calculation is 3,200,000, and the total memory usage is 102,400,000 (Bytes) (3,200,000 ⁇ 4 ⁇ 8: Four variables of 8 bytes per grid (requires x, y, z coordinates and distance function values).
  • the memory usage can be reduced to about 1/40 to 1/100 of the conventional one.
  • the range in which the virtual particles 110 are set (arranged) corresponds to the range (interpolation range) that specifies the virtual particles 110 used to calculate the distance from the particle 101 to the structure 100. It may be a thing. According to this configuration, the required number of virtual particles 110 can be set, and the simulation of the particles 101 and the structure 100 can be performed more efficiently.
  • the structure 100 may receive an interaction force from the particles 101. According to this configuration, since the influence of the particle 101 on the structure 100 can also be simulated, a more accurate simulation of the particle 101 and the structure 100 can be performed.
  • the virtual region may be the virtual particle 110.
  • the virtual area can be handled appropriately and easily.
  • the conventional particle-based simulation method can be easily applied, and the method according to the present embodiment can be speeded up. That is, the simulation of the particle 101 and the structure 100 can be performed appropriately and easily.
  • the virtual region does not necessarily have to be a particulate region, and may have an arbitrary shape.
  • the work area is a three-dimensional space, but may be a space other than three-dimensional, for example, a two-dimensional space.
  • the particle simulation program 30 is inserted into a computer and accessed, or stored in a program storage area 21 formed on a recording medium 20 provided in the computer.
  • the particle simulation program 30 includes a simulation information holding module 31, a structure initial information acquisition module 32, a virtual region setting module 33, a position information acquisition module 34, a virtual region identification module 35, a particle distance calculation module 36, An interaction force calculation module 37 and a time transition module 38 are provided.
  • the function realized by executing the time transition module 38 include the simulation information holding unit 11, the structure initial information acquisition unit 12, the virtual region setting unit 13, and the position information of the particle simulation device 10 described above.
  • the acquisition unit 14, the virtual region specification unit 15, the particle distance calculation unit 16, the interaction force calculation unit 17, and the time transition unit 18 are the same.
  • a part or all of the particle simulation program 30 may be transmitted via a transmission medium such as a communication line and received and recorded (including installation) by another device. Further, each module of the particle simulation program 30 may be installed in any one of a plurality of computers instead of one computer. In that case, the series of particle simulation programs 30 described above is performed by the computer system of the plurality of computers.
  • DESCRIPTION OF SYMBOLS 10 Particle simulation apparatus, 11 ... Simulation information holding part, 12 ... Structure initial information acquisition part, 13 ... Virtual area setting part, 14 ... Position information acquisition part, 15 ... Virtual area specific

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • Mathematical Physics (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Optimization (AREA)
  • Mathematical Analysis (AREA)
  • Pure & Applied Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • Operations Research (AREA)
  • Algebra (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

粒子と構造物との間の距離計算を少ない計算負荷で正確に行い、これにより、粒子と構造物とのシミュレーションを適切に行う。　粒子シミュレーション装置１０は、粒子及び構造物の挙動をシミュレーションする装置である。粒子シミュレーション装置１０は、構造物の初期位置及び形状を示す情報を取得する構造物初期情報取得部１２と、当該情報に基づいて構造物の表面の近傍に複数の仮想粒子を設定し、複数の仮想粒子それぞれについて、構造物までの距離を算出する仮想領域設定部１３と、粒子の位置から予め設定した距離内の仮想粒子を特定する仮想領域特定部１５と、特定された仮想粒子について算出された距離に基づいて、粒子から構造物までの距離を算出する粒子距離算出部１６と、算出された距離に基づいて粒子と構造物との間の相互作用力を算出する相互作用力算出部１７とを備える。

Description

粒子シミュレーション装置、粒子シミュレーション方法及び粒子シミュレーションプログラム
　本発明は、粒子及び構造物の挙動をシミュレーションする粒子シミュレーション装置、粒子シミュレーション方法及び粒子シミュレーションプログラムに関する。
　従来、ＤＥＭ（Discrete Element Method）、ＳＰＨ（Smoothed Particle Hydrodynamics）及びＭＰＳ（Moving Particle Semi-implicit method）等の粒子法シミュレーションによって、粒子と複雑な形状の構造物との挙動の解析が行われている。このようなシミュレーションでは、粒子と構造物との接触判定を行う必要がある。最近では、計算負荷やシミュレーションの正確性等が考慮されて、符号付き距離関数を用いた接触判定が提案されている（例えば、非特許文献１参照）。符号付き距離関数とは、シミュレーションの対象となる空間内を任意サイズの直交格子で分割し、格子点上から構造物表面までの最近接距離を格子点に記憶させたものである（例えば、非特許文献２参照）。粒子から構造物までの距離が、粒子が位置する格子及び隣接格子の距離関数データを用いて内挿補間して求められ、算出された距離に基づいて接触判定が行われる。
Y. Shigeto, M. Sakai, Arbitrary-shaped wall boundary modeling based on signed distance functions for granular flow simulations, Chem. Eng. J. 231 (2013) 464-476. Y. Wu, J. Man, Z. Xie, A double layer method for constructing signed distance fields from triangle meshes, Graphical Models 76 (2014) 214-223.
　構造物の時間変化に伴う移動に対して、距離関数は、例えば、次の何れかのように扱われる。（１）距離関数を固定し、構造物の移動及び回転に対応させて粒子の座標系を平行移動及び回転させる。しかし、この方法は、構造物が複数存在する場合、即ち、移動軸や回転軸が複数存在する場合には適用できない。（２）構造物が動くたびに距離関数を作製し直す。しかし、この方法は距離関数を作り直すための計算コストが膨大となり、実際の計算に用いるのは容易ではない。（３）微分方程式を用いた距離関数の移流計算を行う。しかし、この方法は、格子法であるために数値拡散が生じて距離関数がなまるため、通常、一定時間間隔毎に距離関数を新たに計算し直す再初期化が必要となる。また、その際、それまで格子点上で滑らかに変化していた距離関数が不連続に変化し、不連続に変化する距離振動が発生する。上記の数値拡散や数値振動により、シミュレーションの正確性が低くなり、産業や学術分野では計算精度の面から実用に向かない。（４）構造物の移動に応じて距離関数の時系列データを予め複数用意しておく。しかし、この方法は、格子数×時系列数だけデータ量が必要となり、予め記憶するデータ量が膨大になる。
　本発明は、上記の問題点に鑑みてなされたものであり、粒子と構造物との間の距離計算を少ない計算負荷で正確に行い、これにより、粒子と構造物とのシミュレーションを適切に行うことができる粒子シミュレーション装置、粒子シミュレーション方法及び粒子シミュレーションプログラムを提供することを目的とする。
　上記の目的を達成するために、本発明の一実施形態に係る粒子シミュレーション装置は、作業空間内の粒子と構造物との間の相互作用力に基づいて、当該粒子の位置及び速度を算出し、当該粒子及び当該構造物の挙動をシミュレーションする粒子シミュレーション装置であって、構造物の初期位置及び形状を示す構造物初期情報を取得する構造物初期情報取得手段と、構造物初期情報取得手段によって取得された構造物初期情報によって示される構造物の初期位置及び形状に基づいて、当該構造物の表面の近傍に複数の仮想領域を設定し、設定した複数の仮想領域それぞれについて、構造物までの距離を算出する仮想領域設定手段と、粒子の位置を示す位置情報を取得する位置情報取得手段と、位置情報取得手段によって取得された位置情報によって示される粒子の位置から、予め設定した距離内の仮想領域を特定する仮想領域特定手段と、仮想領域特定手段によって特定された仮想領域について、仮想領域設定手段によって算出された距離に基づいて、粒子から構造物までの距離を算出する粒子距離算出手段と、粒子距離算出手段によって算出された距離に基づいて、粒子と構造物との間の相互作用力を算出する相互作用力算出手段と、相互作用力算出手段によって算出された相互作用力に基づいて、次の時間ステップの粒子の位置及び速度を算出すると共に、当該構造物との位置関係を維持したまま仮想領域を移動させる時間遷移手段と、を備える。
　本発明の一実施形態に係る粒子シミュレーション装置では、構造物の表面の近傍に設定された複数の仮想領域について算出された距離に基づいて、粒子と構造物との間の距離計算が行われる。シミュレーションにおいて、当該仮想領域は、構造物との位置関係を維持したまま移動される。従って、仮想領域についての距離（距離関数）を計算しなおす必要がない。また、従来の距離関数とは異なり構造物の近傍に設定された仮想領域に対してのみ、距離が計算されていればよい。また、仮想領域は、構造物毎に設定されるため、構造物が複数存在する場合、即ち、互いに異なった移動軸や回転軸が複数存在する場合であっても、適用することができる。また、距離に対応付けられた仮想領域が動くため、移流計算による数値拡散も生じない。即ち、本発明の一実施形態に係る粒子シミュレーション装置によれば、粒子と構造物との間の距離計算を少ない計算負荷で正確に行い、これにより、粒子と構造物とのシミュレーションを適切に行うことができる。
　仮想領域設定手段は、予め設定した距離に応じた範囲の構造物の表面の近傍に複数の仮想領域を設定することとしてもよい。この構成によれば、仮想領域を必要な分だけ設定することができ、更に効率的に粒子と構造物とのシミュレーションを行うことができる。
　時間遷移手段は、相互作用力算出手段によって算出された相互作用力に基づいて、次の時間ステップの構造物の位置及び速度を算出することとしてもよい。この構成によれば、構造物に対する粒子の影響もシミュレーションすることができるため、より正確な粒子と構造物とのシミュレーションを行うことができる。
　仮想領域は、粒子状の領域であることとしてもよい。この構成によれば、仮想領域を適切かつ容易に扱うことができ、適切かつ容易に粒子と構造物とのシミュレーションを行うことができる。
　ところで、本発明は、上記のように粒子シミュレーション装置の発明として記述できる他に、以下のように粒子シミュレーション方法及び粒子シミュレーションプログラムの発明としても記述することができる。これはカテゴリが異なるだけで、実質的に同一の発明であり、同様の作用及び効果を奏する。
　即ち、本発明の一実施形態に係る粒子シミュレーション方法は、作業空間内の粒子と構造物との間の相互作用力に基づいて、当該粒子の位置及び速度を算出し、当該粒子及び当該構造物の挙動をシミュレーションする粒子シミュレーション装置の動作方法である粒子シミュレーション方法であって、構造物の初期位置及び形状を示す構造物初期情報を取得する構造物初期情報取得ステップと、構造物初期情報取得ステップにおいて取得された構造物初期情報によって示される構造物の初期位置及び形状に基づいて、当該構造物の表面の近傍に複数の仮想領域を設定し、設定した複数の仮想領域それぞれについて、構造物までの距離を算出する仮想領域設定ステップと、粒子の位置を示す位置情報を取得する位置情報取得ステップと、位置情報取得ステップにおいて取得された位置情報によって示される粒子の位置から、予め設定した距離内の仮想領域を特定する仮想領域特定ステップと、仮想領域特定ステップにおいて特定された仮想領域について、仮想領域設定ステップにおいて算出された距離に基づいて、粒子から構造物までの距離を算出する粒子距離算出ステップと、粒子距離算出手段において算出された距離に基づいて、粒子と構造物との間の相互作用力を算出する相互作用力算出ステップと、相互作用力算出ステップにおいて算出された相互作用力に基づいて、次の時間ステップの粒子の位置及び速度を算出すると共に、当該構造物との位置関係を維持したまま仮想領域を移動させる時間遷移ステップと、を含む。
　また、本発明の一実施形態に係る粒子シミュレーションプログラムは、コンピュータを、作業空間内の粒子と構造物との間の相互作用力に基づいて、当該粒子の位置及び速度を算出し、当該粒子及び当該構造物の挙動をシミュレーションする粒子シミュレーション装置として機能させる粒子シミュレーションプログラムであって、コンピュータを、構造物の初期位置及び形状を示す構造物初期情報を取得する構造物初期情報取得手段と、構造物初期情報取得手段によって取得された構造物初期情報によって示される構造物の初期位置及び形状に基づいて、当該構造物の表面の近傍に複数の仮想領域を設定し、設定した複数の仮想領域それぞれについて、構造物までの距離を算出する仮想領域設定手段と、粒子の位置を示す位置情報を取得する位置情報取得手段と、位置情報取得手段によって取得された位置情報によって示される粒子の位置から、予め設定した距離内の仮想領域を特定する仮想領域特定手段と、仮想領域特定手段によって特定された仮想領域について、仮想領域設定手段によって算出された距離に基づいて、粒子から構造物までの距離を算出する粒子距離算出手段と、粒子距離算出手段によって算出された距離に基づいて、粒子と構造物との間の相互作用力を算出する相互作用力算出手段と、相互作用力算出手段によって算出された相互作用力に基づいて、次の時間ステップの粒子の位置及び速度を算出すると共に、当該構造物との位置関係を維持したまま仮想領域を移動させる時間遷移手段と、として機能させる。
　本発明の一実施形態では、仮想領域についての距離（距離関数）を計算しなおす必要がない。また、従来の距離関数とは異なり構造物の近傍に設定された仮想領域に対してのみ、距離が計算されていればよい。また、仮想領域は、構造物毎に設定されるため、構造物が複数存在する場合、即ち、互いに異なった移動軸や回転軸が複数存在する場合であっても、適用することができる。また、距離に対応付けられた仮想領域が動くため、移流計算による数値拡散も生じない。即ち、本発明の一実施形態によれば、粒子と構造物との間の距離計算を少ない計算負荷で正確に行い、これにより、粒子と構造物とのシミュレーションを適切に行うことができる。
本発明の実施形態に係る粒子シミュレーション装置の機能構成を示す図である。 構造物に対する仮想粒子を模式的に示す図である。 構造物に対する仮想粒子が配置される範囲及び配置された状態を示す図である。 シミュレーション対象の粒子から設定した距離内の仮想粒子を模式的に示す図である。 本発明の実施形態に係る粒子シミュレーション装置で実行される処理（粒子シミュレーション方法）を示すフローチャートである。 本発明の実施形態におけるシミュレーションの実行中の状態（スナップショット）を示す図である。 本発明の実施形態に係る粒子シミュレーションプログラムの構成を、記録媒体と共に示す図である。
　以下、図面と共に本発明に係る粒子シミュレーション装置、粒子シミュレーション方法及び粒子シミュレーションプログラムの実施形態について詳細に説明する。なお、図面の説明においては同一要素には同一符号を付し、重複する説明を省略する。
　図１に本実施形態に係る粒子シミュレーション装置１０を示す。粒子シミュレーション装置１０は、作業空間（計算領域）内の複数の球形（中心から表面までの距離が一定の形状）の粒子及び構造物の挙動をシミュレーション（解析）する装置である。具体的には、粒子シミュレーション装置は、シミュレーション上の時刻である時間ステップ毎に各粒子及び構造物の位置及び速度に基づき、各粒子に働く力を算出する。各粒子及び構造物に働く力には、例えば、粒子と構造物との間、及び粒子間での相互作用である接触（衝突）による相互作用力である接触力が含まれる。なお、粒子と構造物との間、及び粒子間での相互作用は、必ずしも接触に限られず、位置及び速度に影響を及ぼすものであれば任意のものとしてもよい。粒子シミュレーション装置１０は、算出した力に基づいて、次の時間ステップにおける各粒子及び構造物の位置及び速度を算出する。なお、本実施形態における粒子シミュレーションは、ＤＥＭに基づいて行われるが、ＳＰＨやＭＰＳ等の他の粒子法に基づいて行われてもよい。また、粒子法の種類によっては、粒子径や回転速度を情報として保持しなくてもよい。
　本実施形態に係る粒子シミュレーション装置１０によるシミュレーションの対象となる粒子は、従来の粒子シミュレーションの対象となっていた任意の粒子を含む。例えば、土砂や粉体を対象とすることができる。あるいは、流体や固体を複数の粒子からなるものと仮定して対象とすることとしてもよい。本実施形態に係る粒子シミュレーション装置１０によるシミュレーションの対象となる構造物は、任意の物体を含む。例えば、粉体を撹拌、混合する装置（の粉体を撹拌、混合する部分（撹拌翼））を対象とする。シミュレーションの対象となる構造物は、複数あってもよい。
　本実施形態に係る粒子シミュレーション装置１０によるシミュレーションにより、物理的な問題をシミュレーションすることができる。例えば、粉体の撹拌、混合についてのシミュレーションを行うことができる。あるいは、ドリルビットの形状に応じた掘削能力の検証やブルドーザー等の建機や重機の設計にも用いることも考えられる。このように、本実施形態に係る粒子シミュレーション装置１０によるシミュレーションは、土木、粉体等の産業分野での適用が考えられる。また、堤防、建物、地形と津波や土砂との相互作用や土砂崩れ、雪崩等の自然現象のシミュレーションにも用いることができる。このように、本実施形態に係る粒子シミュレーション装置１０によるシミュレーションは、防災分野での適用も考えられる。
　粒子シミュレーション装置１０は、例えば、ＣＰＵ（Central Processing Unit）、ＧＰＵ（Graphics Processing Unit）、メモリ、ハードディスク、ディスプレイ等のハードウェアを備えるコンピュータとして構成される。これらの構成要素がプログラム等により動作することによって、後述する粒子シミュレーション装置１０としての機能が発揮される。粒子シミュレーション装置１０は、並列演算が可能な装置において特に効果的に動作する。なお、粒子シミュレーション装置１０は、演算装置としては必ずしもＧＰＵを備えている必要はなく、ＣＰＵのみを備えた構成であってもよい。
　図１に示すように粒子シミュレーション装置１０は、シミュレーション情報保持部１１と、構造物初期情報取得部１２と、仮想領域設定部１３と、位置情報取得部１４と、仮想領域特定部１５と、粒子距離算出部１６と、相互作用力算出部１７と、時間遷移部１８とを備えて構成される。
　本実施形態において粒子が運動する領域である作業領域は、例えば、三次元の空間であり、一辺の大きさが予め設定された立方体のセル（メッシュ）に分割（区分）されている。なお、セルは立方体以外の形状、例えば、直方体であってもよい。粒子シミュレーション装置１０は、シミュレーションの処理を行う前に作業空間を予めセルに分割しており、作業空間がどのようにセルに分割されているか予め把握している。上記の一辺の大きさは、例えば、シミュレーションの対象となる複数の粒子の粒子径等に応じて予め設定されている。また、作業領域内の各セルにはセルを特定するセル番号が付されている。セル番号は、例えば、作業空間内のセルの位置に応じて順番に付されている。
　シミュレーション情報保持部１１は、シミュレーションに用いられる情報を保持する手段である。具体的には、シミュレーション情報保持部１１は、作業領域内の複数の粒子それぞれについての粒子情報を保持する。粒子情報は、粒子の座標、粒子の速度及び粒子半径を示す情報を含む。粒子の座標は、作業空間における粒子の位置を示す三次元座標である。粒子の速度は、並進速度及び回転速度を含む。粒子の座標及び粒子の速度については、シミュレーションの開始時の情報（初期情報）は、予めシミュレーション情報保持部１１に粒子シミュレーション装置１０のユーザ等により入力されており、また、シミュレーション中の情報は、後述する時間遷移部１８によって更新される。粒子半径は、予めシミュレーション情報保持部１１に粒子シミュレーション装置１０のユーザ等により入力されている。
　シミュレーション情報保持部１１は、作業領域内の構造物それぞれについての構造物情報を保持する。構造物情報は、構造物の座標、構造物の速度及び構造物の形状を示す情報を含む。構造物の座標は、作業空間における構造物の位置（例えば、構造物の特定の箇所の位置）を示す三次元座標である。構造物の速度は、並進速度及び回転速度を含む。構造物の座標及び構造物の速度については、粒子の座標及び粒子の速度と同様に、シミュレーションの開始時の情報（初期情報）は、予めシミュレーション情報保持部１１に粒子シミュレーション装置１０のユーザ等により入力されており、また、シミュレーション中の情報は、後述する時間遷移部１８によって更新される。構造物の形状を示す情報は、例えば、ポリゴンで構造物の形状を表したＣＡＤ（computer-aided design）データである。構造物の形状を示す情報は、予めシミュレーション情報保持部１１に粒子シミュレーション装置１０のユーザ等により入力されている。構造物情報は、構造物がどのように移動するか（例えば、移動軸や回転軸等）を示す情報を含んでいてもよい。当該情報も、予めシミュレーション情報保持部１１に粒子シミュレーション装置１０のユーザ等により入力されている。
　また、シミュレーション情報保持部１１は、後述するように仮想粒子の情報を保持する。また、シミュレーション情報保持部１１は、上記の粒子情報及び構造物情報以外のシミュレーションに利用される情報を、予め入力して保持していてもよい。このような情報としては、摩擦係数、弾性係数、粘性減衰係数、反発係数等である。
　構造物初期情報取得部１２は、構造物の初期位置及び形状を示す構造物初期情報を取得する構造物初期情報取得手段である。例えば、構造物初期情報取得部１２は、シミュレーション情報保持部１１に保持されている構造物に係る初期情報、及び構造物の形状を示す情報を構造物初期情報として取得する。構造物初期情報取得部１２は、取得した構造物初期情報を仮想領域設定部１３に出力する。
　仮想領域設定部１３は、構造物初期情報取得部１２によって取得された構造物初期情報によって示される構造物の初期位置及び形状に基づいて、当該構造物の表面の近傍に複数の仮想領域を設定し、設定した複数の仮想領域それぞれについて、構造物までの距離（最近接距離）を算出する仮想領域設定手段である。仮想領域は、後述するように、シミュレーション対象の粒子から構造物までの距離の算出に用いる領域である。本実施形態では、仮想領域は、球形（中心から表面までの距離が一定の形状）の粒子状の領域であり、仮想粒子と呼ぶ。複数の仮想粒子の径の長さは一律であり、予め設定されている。仮想粒子は、計算の負荷等を考慮して、シミュレーション対象の粒子よりも大きいものとしてもよい。
　上述したように、従来は、シミュレーション対象の粒子から構造物までの距離は、非特許文献１に示されるように符号付き距離関数が用いられて求められていた。例えば、符号付き距離関数は、作業領域がセルに（直交格子で）分割された格子点に構造物からの距離を記憶させたものである。図２（ａ）に格子点状の距離関数値を基にして生成された構造物１００からの等しい距離の位置を線で結んだ等値線図を示す。本実施形態では、図２（ｂ）に示すように、構造物１００の表面（構造物１００とそれ以外の部分との境界）の近傍に複数の仮想粒子１１０を設定する。
　具体的には、仮想領域設定部１３は、まず、構造物初期情報によって初期位置及び形状が示される構造物１００についての符号付き距離関数を算出する。この算出は、非特許文献２に記載されているような従来の方法を用いて行えばよい。仮想領域設定部１３は、算出した符号付き距離関数に基づいて、作業領域内の構造物１００からの近傍の位置を検出する。構造物からの近傍は、図２（ｂ）に示すように、構造物の外側、表面及び内側の領域を含む。仮想領域設定部１３は、当該近傍の領域に、互いに接触するように複数の仮想粒子１１０を設定する。例えば、作業領域が三次元である場合には、仮想領域設定部１３は、当該領域内において単純立方格子又は面心立方格子の格子点を設定し、当該格子点に仮想粒子１１０を設定（配置）する。なお、仮想粒子の設定には、構造物１００からの距離が分かればよいので必ずしも符号付き距離関数を用いる必要はなく、任意の方法で行われてもよい。
　仮想粒子１１０を設定する構造物１００の表面の部分は、シミュレーション対象の粒子が接触する（相互作用する）可能性がある部分である。構造物１００の表面全体に渡って粒子が接触する可能性がある場合には、図２（ｂ）に示すように、構造物１００の表面全体に渡って表面の近傍に仮想粒子１１０を設定する。接触する可能性がない場所が予め分かっている場合には、その部分の近傍には仮想粒子１１０を設定しなくてもよい。
　仮想領域設定部１３は、構造物１００の表面の法線方向においては、所定の距離（所定の厚み）の範囲を構造物からの近傍とする。例えば、図３に示すように、構造物１００の表面から外側へ距離ｄ離れた位置をその距離の中心として、そこから外側及び内側（図３における上下方向）へそれぞれ距離Ｈ分の領域を構造物１００の近傍の領域とする。距離ｄは、構造物１００と、シミュレーション対象の粒子１０１とが相互作用し始める、構造物１００の表面と粒子１０１の中心との間の距離である。なお、シミュレーション対象の粒子１０１の大きさが一律でない等、相互作用し始める距離が粒子１０１毎に一律でない場合には、例えば、図３におけるｄ－Ｈについてはそれらの最小値、ｄ＋Ｈについてはそれらの最大値を用いることとしてもよい。
　相互作用が接触である場合、距離ｄは、粒子１０１の半径である。距離Ｈは、粒子１０１から構造物１００までの距離を算出するために用いる仮想粒子１１０を特定するための距離である（詳細には後述する）。距離Ｈは、予め設定されている。上記の通り、仮想領域設定部１３は、距離Ｈに応じた範囲の構造物１００の表面の近傍に複数の仮想粒子１１０を設定することとしてもよい。なお、距離ｄが距離Ｈに比べて大きい場合、構造物１００の表面及び内側には、仮想粒子１１０が設定されない場合もある。また、仮想粒子１１０が設定される領域（法線方向の距離、厚み）は、上記のように設定（算出）される必要はなく、予め設定されていてもよい。
　後述するように、構造物１００と、設定された仮想粒子１１０とは、シミュレーションの間中、位置関係が固定される。
　仮想領域設定部１３は、設定した各仮想粒子１１０について、当該仮想粒子１１０から（配置された表面に係る）構造物１００までの距離、例えば、仮想粒子１１０の中心の位置から構造物１００の表面までの距離を算出する。この算出には、従来の符号付き距離関数を用いた内挿補間による方法等の任意の方法を用いることができる。なお、仮想粒子１１０から構造物１００までの距離の算出は、仮想粒子１１０の配置とあわせて行われてもよい。設定された仮想粒子１１０は、仮想粒子１１０毎に一意となる番号を付与する等して、一意に特定できるようにしておく。仮想領域設定部１３は、設定した仮想粒子１１０の位置及び上記の距離を示す情報をシミュレーション情報保持部１１に格納して保持させる。
　仮想領域設定部１３は、仮想粒子１１０の距離の算出の際に、構造物１００の表面に対する法線ベクトルを算出して、距離とあわせて対応付けておいてもよい。法線ベクトルの情報は、後述する内挿補間の処理において用いられ得る。また、仮想粒子１１０には、ヤング率等の機械的物性を対応付けておき、後述する処理に用いることもできる。
　上記の構造物初期情報取得部１２及び仮想領域設定部１３による処理は、シミュレーションの各時間ステップに係る演算の前に行われる。以降の機能部は、シミュレーションの各時間ステップに係る演算のための構成である。
　位置情報取得部１４は、シミュレーション対象となる粒子１０１の位置を示す位置情報を取得する位置情報取得手段である。具体的には、位置情報取得部１４は、シミュレーション情報保持部１１に保持された各粒子１０１の現時間ステップの粒子情報を、粒子１０１の位置情報として取得する。位置情報取得部１４は、取得した粒子１０１の粒子情報を仮想領域特定部１５に入力する。
　位置情報取得部１４は、仮想粒子１１０の位置を示す位置情報も取得する。位置情報取得部１４は、シミュレーション情報保持部１１に保持された各粒子１０１の現時間ステップの位置を示す情報を取得する。位置情報取得部１４は、取得した仮想粒子１１０の位置情報を仮想領域特定部１５に出力する。
　仮想領域特定部１５は、位置情報取得部１４によって取得された位置情報によって示される粒子１０１の位置から、予め設定した距離Ｈ内の仮想粒子１１０を特定する仮想領域特定手段である。ここで特定された仮想粒子１１０は、粒子１０１から構造物１００までの距離を内挿補間によって算出するために用いられる。即ち、距離Ｈは、内挿補間を行う範囲（補間範囲）を示す距離であり、内挿補間によって適切に粒子１０１から構造物１００までの距離を算出できるような値が予め設定される。具体的には、距離Ｈは、仮想粒子１１０が単純立方格子状に設定される場合には、Ｈ＞Ｄ√ｎを満たすように、仮想粒子１１０が面心立方格子状に設定される場合には、Ｈ＞Ｄを満たすようにされるのがよい。なお、上記の関係式で、Ｄは図３に示すように仮想粒子１１０の直径であり、ｎは作業領域の次元数である。
　図４に示すように、仮想領域特定部１５は、粒子１０１毎に、位置情報取得部１４から入力した情報によって示される粒子１０１（粒子ｉ）の位置（粒子１０１の中心の位置）と各仮想粒子１１０（粒子ｊ）の位置（仮想粒子１１０の中心の位置）との間の距離Ｐｉｊを算出する。仮想領域特定部１５は、算出した距離Ｐｉｊが距離Ｈ以内であるか否かを判断する。仮想領域特定部１５は、当該判断に基づき、粒子１０１毎に、粒子１０１からの距離が距離Ｈ以内である仮想粒子１１０を特定する。図４に示す例では、４つの仮想粒子１１０ａが特定される。仮想領域特定部１５は、粒子１０１毎に特定された仮想粒子１１０を示す情報を粒子距離算出部１６に出力する。
　なお、距離Ｈ内の仮想粒子１１０の特定は、従来の粒子シミュレーションにおける粒子間の接触判定（粒子間距離が閾値以内かの判定）の方法を利用して、効率的に行うこととしてもよい。例えば、特開２０１５－１１５５６７号公報（特許文献１）及び特開２０１０－２３８０３０号公報（特許文献２）に記載されたような、粒子１０１及び仮想粒子１１０が位置するセルに基づいて、粒子１０１と仮想粒子１１０とのペアを設定し、当該ペアについて距離の判断を行う方法をとってもよい。
　粒子距離算出部１６は、仮想領域特定部１５によって特定された仮想粒子１１０について、仮想領域設定部１３によって算出された距離に基づいて、粒子１０１から構造物１００までの距離（最近接距離）を算出する粒子距離算出手段である。粒子距離算出部１６は、シミュレーション情報保持部１１に保持された、粒子１０１毎に特定された仮想粒子１１０についての構造物１００までの距離を示す情報を取得する。粒子距離算出部１６は、取得した情報によって示される距離に基づいて、粒子１０１から構造物１００までの距離（図４に示す粒子１０１（粒子ｉ）の中心の位置から構造物１００の表面までの距離ｄｉｓｉを算出する。この算出は、例えば、非特許文献１に示される従来の内挿補間によって行われてもよいし、その他の任意の方法で行うこととしてもよい。例えば、ＳＰＨ等で用いられるカーネル近似も用いることができる。
　特定された仮想粒子１１０がない、あるいは、特定された仮想粒子１１０が内挿補間を行うのに不十分な数である（例えば、予め設定した閾値以下である）場合には、当該粒子１０１については、構造物１００までの距離が少なくとも相互作用が発生する距離にはなっていないとして、構造物１００までの距離を算出しない。この場合、当該距離を用いた処理も行われない。粒子距離算出部１６は、粒子１０１毎に算出した構造物１００までの距離を示す情報を相互作用力算出部１７に出力する。
　また、粒子距離算出部１６は、通常、上記の距離算出の際に粒子１０１から構造物１００の表面に対する法線ベクトルを算出する。法線ベクトルによって、相互作用力算出部１７によって算出される相互作用力をベクトル量で表すことができる。当該ベクトル量が表す方向は、相互作用力が作用する方向を表す。法線ベクトルは、粒子１０１毎に特定された仮想粒子１１０についての構造物１００までの距離（距離関数データ）に対して空間勾配を計算することで求めることができる。なお、空間勾配からの法線ベクトルの算出は、例えば、非特許文献１に示される従来の方法によって行われてもよい。あるいは、仮想粒子１１０について予め法線ベクトルを算出している場合には、その情報を内挿補間して算出することとしてもよい。その他の任意の方法で行うこととしてもよい。例えば、ＳＰＨ等で用いられるカーネル近似による勾配算出法も用いることができる。粒子距離算出部１６は、算出した法線ベクトルを示す情報を、距離を示す情報とあわせて相互作用力算出部１７に出力する。
　相互作用力算出部１７は、粒子距離算出部１６によって算出された距離に基づいて、粒子１０１と構造物１００との間の相互作用力を算出する相互作用力算出手段である。相互作用力算出部１７は、粒子１０１から構造物１００までの距離が、予め設定された閾値よりも小さいか否かを判断することで粒子１０１と構造物１００との間で相互作用が発生し得るか否かの判定（相互作用判定）を行う。相互作用力算出部１７は、当該距離が閾値よりも小さいと判断した場合、粒子１０１と構造物１００との間で相互作用が発生すると判定する。相互作用力算出部１７は、当該距離が閾値よりも小さくないと判断した場合、粒子１０１と構造物１００との間で相互作用が発生しないと判定する。
　例えば、相互作用が接触である場合、相互作用力算出部１７は、当該距離が閾値である粒子１０１の半径ｄよりも小さいか否かを判断し、当該距離が閾値よりも小さいと判断した場合、粒子１０１と構造物１００とが接触しており、それらの間に接触力が発生すると判定する。相互作用力の計算は、従来の粒子シミュレーションと同様の方法やその他の任意の方法で行うことができる。相互作用力算出部１７は、通常、粒子距離算出部１６によって算出された粒子１０１から構造物１００の表面に対する法線ベクトルに基づいて、相互作用力をベクトル量として算出する。
　相互作用力算出部１７は、粒子１０１と構造物１００との間の相互作用力だけでなく、粒子１０１同士の相互作用力、及び構造物１００同士の相互作用力を算出してもよい。粒子１０１同士の相互作用力の計算は、特許文献１及び２に記載された従来の粒子シミュレーションと同様の方法やその他の任意の方法で行うことができる。構造物１００同士の相互作用力の計算も、従来のシミュレーションと同様の方法やその他の任意の方法で行うことができる。なお、相互作用力の算出には、シミュレーション情報保持部１１に保持されている、相互作用力の計算に必要なパラメータを取得し、相互作用力の計算に用いてもよい。
　相互作用力算出部１７は、上記のように算出した個々の粒子１０１と構造物１００との間の相互作用力、粒子１０１同士の相互作用力、及び構造物１００同士の相互作用力から、各粒子１０１及び各構造物１００についての相互作用力の総和を算出する。相互作用力算出部１７は、算出した各粒子１０１及び各構造物１００の相互作用力の総和を示す情報を時間遷移部１８に出力する。
　時間遷移部１８は、相互作用力算出部１７によって算出された相互作用力に基づいて、次の時間ステップの粒子１０１及び構造物１００の位置及び速度を算出すると共に、当該構造物１００との位置関係を維持したまま仮想粒子１１０を移動させる時間遷移手段である。具体的には、時間遷移部１８は、粒子１０１及び構造物１００毎にシミュレーション情報保持部１１に保持された粒子情報及び構造物情報を取得し、粒子情報及び構造物情報に示される現時間ステップの座標及び速度と相互作用力の総和とから、次の時間ステップの座標及び速度を算出する。この算出は、例えば、従来のシミュレーションと同様の方法やその他の任意の方法で行うことができる。時間遷移部１８は、算出した次の時間ステップにおける粒子１０１及び構造物１００の位置及び速度で、粒子１０１及び構造物１００毎にシミュレーション情報保持部１１に保持された粒子情報及び構造物情報を更新する。
　時間遷移部１８は、構造物１００について次の時間ステップの位置を算出すると、当該構造物１００に対して設定された仮想粒子１１０を、当該構造物１００との位置関係を維持したまま（当該構造物１００の移動と同様に）移動させる。時間遷移部１８は、移動させた仮想粒子１１０の位置を示す情報で、シミュレーション情報保持部１１に保持された仮想粒子１１０を示す情報を更新する（当該情報を、次の時間ステップの仮想粒子１１０の位置を示す情報として、シミュレーション情報保持部１１に保持させる）。
　なお、構造物１００については、相互作用力に基づいた位置及び速度の算出を行わなくてもよい。例えば、構造物１００に対する粒子１０１の影響が極めて小さく、また、構造物１００同士の相互作用が生じない場合については、当該算出を行わなくてもよい。この場合、例えば、構造物１００が予めどのように動くかを示す情報が、シミュレーション情報保持部１１に記憶されており、時間遷移部１８は、それに応じて構造物１００の移動（次の時間ステップの位置及び速度）を算出することとしてもよい。また、時間遷移部１８は、構造物１００の当該移動に応じて（予め決まっている構造物１００の表面の軌跡に応じて）、仮想粒子１１０を移動させる。
　時間遷移部１８によって、全ての粒子１０１、構造物１００及び仮想粒子１１０についての更新が行われると、次の時間ステップの処理が行われる。
　粒子シミュレーション装置１０では、１つの時間ステップでの演算が完了する毎に、シミュレーションの終了条件を満たしているか否かが判断される。例えば、予め設定した回数（時間ステップ）の演算が終了した場合、終了条件を満たしていると判断される。終了条件を満たしていると判断された場合には、粒子シミュレーション装置１０では、シミュレーションが終了される。この場合、例えば、表示装置や他の装置への演算結果の出力等が行われる。終了条件を満たしていないと判断された場合には、次の時間ステップの演算が繰り返し行われる。以上が、粒子シミュレーション装置１０の構成である。
　引き続いて、図５のフローチャートを用いて、本実施形態に係る粒子シミュレーション装置１０の動作方法である、粒子シミュレーション装置１０で実行される処理（粒子シミュレーション方法）を説明する。本処理は、例えば、粒子シミュレーション装置１０のユーザが、粒子シミュレーション装置１０に対してシミュレーションを開始する操作を行うことで開始される。
　粒子シミュレーション装置１０では、構造物初期情報取得部１２によって、シミュレーション情報保持部１１に保持されている構造物１００の初期位置及び形状を示す構造物初期情報が取得される（Ｓ０１、構造物初期情報取得ステップ）。構造物１００の形状の例を図６（ａ）に示す。取得された構造物初期情報は、構造物初期情報取得部１２から仮想領域設定部１３に出力される。
　続いて、仮想領域設定部１３によって、構造物初期情報に係る構造物１００の符号付き距離関数が算出される（Ｓ０２、仮想領域設定ステップ）。続いて、仮想領域設定部１３によって、当該符号付き距離関数が用いられて、複数の仮想粒子１１０が構造物１００の表面の近傍に配置（設定）される（Ｓ０３、仮想領域設定ステップ）。構造物１００の表面の近傍に配置された仮想粒子１１０の例を図６（ｂ）に示す。続いて、仮想領域設定部１３によって、仮想粒子１１０について、当該仮想粒子１１０から構造物１００までの距離が算出される（Ｓ０４、仮想領域設定ステップ）。仮想粒子１１０の位置及び距離を示す情報は、シミュレーション情報保持部１１に保持される。以上がシミュレーションの各時間ステップに係る演算の前に行われる処理である。
　引き続いて、位置情報取得部１４によって、シミュレーション情報保持部１１に保持されている、シミュレーション対象となる粒子１０１及び仮想粒子１１０の位置を示す位置情報が取得される（Ｓ０５、位置情報取得ステップ）。取得された位置情報は、位置情報取得部１４から仮想領域特定部１５に出力される。
　続いて、仮想領域特定部１５によって、位置情報によって示される粒子１０１の位置から距離Ｈ内の仮想粒子１１０が特定される（Ｓ０６、仮想領域特定ステップ）。粒子１０１毎に特定された仮想粒子１１０を示す情報は、仮想領域特定部１５から粒子距離算出部１６に出力される。
　続いて、粒子距離算出部１６によって、粒子１０１毎に特定された仮想粒子１１０についての構造物１００までの距離に基づいて、粒子１０１から構造物１００までの距離が算出される（Ｓ０７、粒子距離算出ステップ）。この算出は、上述したように例えば、内挿補間によって行われる。粒子１０１毎に算出された構造物１００までの距離を示す情報は、粒子距離算出部１６から相互作用力算出部１７に出力される。
　続いて、相互作用力算出部１７によって、粒子１０１毎に算出された構造物１００までの距離に基づいて、粒子１０１と構造物１００との間の相互作用力が算出される（Ｓ０８、相互作用力算出ステップ）。この際、粒子１０１同士の相互作用力、及び構造物１００同士の相互作用力が算出されてもよい。続いて、相互作用力算出部１７によって、算出された相互作用力から、各粒子１０１及び各構造物１００についての相互作用力の総和が算出される。算出された各粒子１０１及び各構造物１００の相互作用力の総和を示す情報は、相互作用力算出部１７から時間遷移部１８に出力される。
　続いて、時間遷移部１８によって、算出された相互作用力に基づいて、次の時間ステップの粒子１０１及び構造物１００の位置及び速度が算出されると共に、当該構造物１００との位置関係を維持したまま仮想粒子１１０が移動される（Ｓ０９、時間遷移ステップ）。上記の演算に基づき、シミュレーション情報保持部１１に保持された粒子情報、構造物情報、仮想粒子１１０の位置を示す情報が次の時間ステップの情報に更新される。
　続いて、粒子シミュレーション装置１０では、シミュレーションの終了条件を満たしているか否かが判断される（Ｓ１０）。終了条件を満たしていると判断された場合（Ｓ１０のＹＥＳ）には、処理（シミュレーション）が終了される。終了条件を満たしていないと判断された場合（Ｓ１０のＮＯ）には、時間ステップが一つ進められて、次の時間ステップでの上述した処理（Ｓ０５～Ｓ１０）が行われる。
　図６（ｃ）～図６（ｆ）に各時間ステップにおける状態の例を示す。図６（ｃ）から図６（ｆ）に進むに従って、時間ステップが進んでいる。このシミュレーションでは、容器の内側に粉体を撹拌するための回転する構造物（撹拌翼）が２つ設けられており、粉体が上部から容器の中に落ちて構造物で撹拌される。なお、それぞれの構造物には、それぞれ別の撹拌速度を与えることができる。以上が、本実施形態に係る粒子シミュレーション装置１０で実行される処理である。
　上述したように本実施形態によれば、構造物１００の表面の近傍に設定された複数の仮想粒子１１０について算出された距離に基づいて、粒子１０１と構造物１００との間の距離計算が行われる。シミュレーションにおいて、当該仮想粒子１１０は、構造物１００との位置関係を維持したまま移動される。従って、仮想粒子１１０についての距離（距離関数）を計算しなおす必要がない。
　また、作業空間内全域に渡る格子点での従来の距離関数とは異なり、構造物１００の近傍に設定された仮想粒子１１０に対してのみ、距離が計算されていればよい。また、個々の仮想粒子１１０は、構造物１００毎に設定され、構造物１００毎の動きに伴って独立に動かすことができる。即ち、構造物１００が複数存在する場合、即ち、互いに異なった移動軸や回転軸が複数存在する場合であっても、適用することができる。
　また、距離に対応付けられた仮想粒子１１０が動くため、移流計算による数値拡散も生じず、距離関数の初期化処理は必要ない。また、仮想粒子１１０を動かすだけであるので、移流計算に比べて非常に簡便で計算負荷も少ない。即ち、本実施形態によれば、粒子１０１と構造物１００との間の距離計算を少ない計算負荷で正確に行い、これにより、粒子１０１と構造物１００とのシミュレーションを適切に行うことができる。
　また、上記の通り、構造物１００の近傍に設定された仮想粒子１１０に対してのみ、距離が計算されて記憶されていればよいため、従来の距離関数と比べて記憶すべきデータ量（メモリ使用量）も圧倒的に削減することができる。
　メモリ使用量の例を以下に示す。作業空間のサイズを１６０Ｄ×２００Ｄ×１００Ｄ、補間範囲Ｈを３．５ｄ、仮想粒子１１０の直径（格子サイズ）をＤ、相互作用距離ｄを０．５Ｄとした場合の例を示す。本実施形態による方法では、距離計算に必要な粒子数は８５，３３６であり、総メモリ使用量は２，７３０，７５２（Ｂｙｔｅ）（８５，３３６×４×８：一粒子当り８ｂｙｔｅの変数が４個（ｘ，ｙ，ｚ座標と距離関数値）必要）である。従来の方法（粒子の座標系を動かす場合）では、距離計算に必要な格子数は３，２００，０００であり、総メモリ使用量は１０２，４００，０００（Ｂｙｔｅ）（３，２００，０００×４×８：一格子当り８ｂｙｔｅの変数が４個（ｘ，ｙ，ｚ座標と距離関数値）必要）である。更に、距離関数の移流計算を行う場合、例えば、ＣＩＰ法（Constrained Interpolation Profile Scheme）を用いるとすると、メモリ使用量は１０２，４００，０００＋３，２００，０００×７×８＝２８１，６００，０００（Ｂｙｔｅ）（一格子当り８ｂｙｔｅの変数が７個必要）である。また、距離関数の時系列データを使う場合、メモリ使用量は１０２，４００，０００＋３，２００，０００×８＝１２８，０００，０００（Ｂｙｔｅ）（時間補間を行うため、２時刻の距離関数データが必要で、距離関数値が一格子当り８Ｂｙｔｅであり、全部で３，２００，０００×８Ｂｙｔｅ追加で必要）である。
　上記の例のように本実施形態では、メモリ使用量を従来の約１／４０～１／１００に削減することができる。
　また、本実施形態のように、仮想粒子１１０を設定（配置）する範囲を、粒子１０１から構造物１００までの距離を算出するのに用いる仮想粒子１１０を特定する範囲（補間範囲）に応じたものとしてもよい。この構成によれば、仮想粒子１１０を必要な分だけ設定することができ、更に効率的に粒子１０１と構造物１００とのシミュレーションを行うことができる。
　また、本実施形態のように、構造物１００が、粒子１０１からの相互作用力を受けることとしてもよい。この構成によれば、構造物１００に対する粒子１０１の影響もシミュレーションすることができるため、より正確な粒子１０１と構造物１００とのシミュレーションを行うことができる。
　また、本実施形態のように、仮想領域を仮想粒子１１０とすることとしてもよい。この構成によれば、仮想領域を適切かつ容易に扱うことができる。例えば、上述したように、従来の粒子ベースのシミュレーションの手法を容易に適用することができ、本実施形態による手法を高速化することができる。即ち、適切かつ容易に粒子１０１と構造物１００とのシミュレーションを行うことができる。但し、必ずしも仮想領域を粒子状の領域とする必要はなく、任意の形状の領域としてもよい。
　なお、上述した実施形態では、作業領域は三次元の空間としたが、三次元以外の空間、例えば、二次元の空間であってもよい。
　引き続いて、上述した一連の粒子シミュレーション装置１０による処理をコンピュータに実行させるための粒子シミュレーションプログラムを説明する。図７に示すように、粒子シミュレーションプログラム３０は、コンピュータに挿入されてアクセスされる、あるいはコンピュータが備える記録媒体２０に形成されたプログラム格納領域２１内に格納される。
　粒子シミュレーションプログラム３０は、シミュレーション情報保持モジュール３１と、構造物初期情報取得モジュール３２と、仮想領域設定モジュール３３と、位置情報取得モジュール３４と、仮想領域特定モジュール３５と、粒子距離算出モジュール３６と、相互作用力算出モジュール３７と、時間遷移モジュール３８とを備えて構成される。シミュレーション情報保持モジュール３１と、構造物初期情報取得モジュール３２と、仮想領域設定モジュール３３と、位置情報取得モジュール３４と、仮想領域特定モジュール３５と、粒子距離算出モジュール３６と、相互作用力算出モジュール３７と、時間遷移モジュール３８とを実行させることにより実現される機能は、上述した粒子シミュレーション装置１０のシミュレーション情報保持部１１と、構造物初期情報取得部１２と、仮想領域設定部１３と、位置情報取得部１４と、仮想領域特定部１５と、粒子距離算出部１６と、相互作用力算出部１７と、時間遷移部１８とそれぞれ同様である。
　なお、粒子シミュレーションプログラム３０は、その一部又は全部が、通信回線等の伝送媒体を介して伝送され、他の機器により受信されて記録（インストールを含む）される構成としてもよい。また、粒子シミュレーションプログラム３０の各モジュールは、１つのコンピュータでなく、複数のコンピュータのいずれかにインストールされてもよい。その場合、当該複数のコンピュータによるコンピュータシステムよって上述した一連の粒子シミュレーションプログラム３０の処理が行われる。
　１０…粒子シミュレーション装置、１１…シミュレーション情報保持部、１２…構造物初期情報取得部、１３…仮想領域設定部、１４…位置情報取得部、１５…仮想領域特定部、１６…粒子距離算出部、１７…相互作用力算出部、１８…時間遷移部、２０…記録媒体、２１…プログラム格納領域、３０…粒子シミュレーションプログラム、３１…シミュレーション情報保持モジュール、３２…構造物初期情報取得モジュール、３３…仮想領域設定モジュール、３４…位置情報取得モジュール、３５…仮想領域特定モジュール、３６…粒子距離算出モジュール、３７…相互作用力算出モジュール、３８…時間遷移モジュール。

Claims (6)

  1. 　作業空間内の粒子と構造物との間の相互作用力に基づいて、当該粒子の位置及び速度を算出し、当該粒子及び当該構造物の挙動をシミュレーションする粒子シミュレーション装置であって、
    　構造物の初期位置及び形状を示す構造物初期情報を取得する構造物初期情報取得手段と、
    　前記構造物初期情報取得手段によって取得された構造物初期情報によって示される構造物の初期位置及び形状に基づいて、当該構造物の表面の近傍に複数の仮想領域を設定し、設定した複数の仮想領域それぞれについて、構造物までの距離を算出する仮想領域設定手段と、
    　粒子の位置を示す位置情報を取得する位置情報取得手段と、
    　前記位置情報取得手段によって取得された位置情報によって示される粒子の位置から、予め設定した距離内の仮想領域を特定する仮想領域特定手段と、
    　前記仮想領域特定手段によって特定された仮想領域について、前記仮想領域設定手段によって算出された距離に基づいて、粒子から構造物までの距離を算出する粒子距離算出手段と、
    　前記粒子距離算出手段によって算出された距離に基づいて、粒子と構造物との間の相互作用力を算出する相互作用力算出手段と、
    　前記相互作用力算出手段によって算出された相互作用力に基づいて、次の時間ステップの粒子の位置及び速度を算出すると共に、当該構造物との位置関係を維持したまま前記仮想領域を移動させる時間遷移手段と、
    を備える粒子シミュレーション装置。
  2. 　前記仮想領域設定手段は、前記予め設定した距離に応じた範囲の前記構造物の表面の近傍に複数の仮想領域を設定する請求項１に記載の粒子シミュレーション装置。
  3. 　前記時間遷移手段は、前記相互作用力算出手段によって算出された相互作用力に基づいて、次の時間ステップの構造物の位置及び速度を算出する請求項１又は２に記載の粒子シミュレーション装置。
  4. 　前記仮想領域は、粒子状の領域である請求項１又は２に記載の粒子シミュレーション装置。
  5. 　作業空間内の粒子と構造物との間の相互作用力に基づいて、当該粒子の位置及び速度を算出し、当該粒子及び当該構造物の挙動をシミュレーションする粒子シミュレーション装置の動作方法である粒子シミュレーション方法であって、
    　構造物の初期位置及び形状を示す構造物初期情報を取得する構造物初期情報取得ステップと、
    　前記構造物初期情報取得ステップにおいて取得された構造物初期情報によって示される構造物の初期位置及び形状に基づいて、当該構造物の表面の近傍に複数の仮想領域を設定し、設定した複数の仮想領域それぞれについて、構造物までの距離を算出する仮想領域設定ステップと、
    　粒子の位置を示す位置情報を取得する位置情報取得ステップと、
    　前記位置情報取得ステップにおいて取得された位置情報によって示される粒子の位置から、予め設定した距離内の仮想領域を特定する仮想領域特定ステップと、
    　前記仮想領域特定ステップにおいて特定された仮想領域について、前記仮想領域設定ステップにおいて算出された距離に基づいて、粒子から構造物までの距離を算出する粒子距離算出ステップと、
    　前記粒子距離算出ステップにおいて算出された距離に基づいて、粒子と構造物との間の相互作用力を算出する相互作用力算出ステップと、
    　前記相互作用力算出ステップにおいて算出された相互作用力に基づいて、次の時間ステップの粒子の位置及び速度を算出すると共に、当該構造物との位置関係を維持したまま前記仮想領域を移動させる時間遷移ステップと、
    を含む粒子シミュレーション方法。
  6. 　コンピュータを、作業空間内の粒子と構造物との間の相互作用力に基づいて、当該粒子の位置及び速度を算出し、当該粒子及び当該構造物の挙動をシミュレーションする粒子シミュレーション装置として機能させる粒子シミュレーションプログラムであって、
    　前記コンピュータを、
    　構造物の初期位置及び形状を示す構造物初期情報を取得する構造物初期情報取得手段と、
    　前記構造物初期情報取得手段によって取得された構造物初期情報によって示される構造物の初期位置及び形状に基づいて、当該構造物の表面の近傍に複数の仮想領域を設定し、設定した複数の仮想領域それぞれについて、構造物までの距離を算出する仮想領域設定手段と、
    　粒子の位置を示す位置情報を取得する位置情報取得手段と、
    　前記位置情報取得手段によって取得された位置情報によって示される粒子の位置から、予め設定した距離内の仮想領域を特定する仮想領域特定手段と、
    　前記仮想領域特定手段によって特定された仮想領域について、前記仮想領域設定手段によって算出された距離に基づいて、粒子から構造物までの距離を算出する粒子距離算出手段と、
    　前記粒子距離算出手段によって算出された距離に基づいて、粒子と構造物との間の相互作用力を算出する相互作用力算出手段と、
    　前記相互作用力算出手段によって算出された相互作用力に基づいて、次の時間ステップの粒子の位置及び速度を算出すると共に、当該構造物との位置関係を維持したまま前記仮想領域を移動させる時間遷移手段と、
    として機能させる粒子シミュレーションプログラム。
PCT/JP2017/008168 2016-03-03 2017-03-01 粒子シミュレーション装置、粒子シミュレーション方法及び粒子シミュレーションプログラム WO2017150626A1 (ja)

Priority Applications (4)

Application Number Priority Date Filing Date Title
CN201780014249.XA CN108701275A (zh) 2016-03-03 2017-03-01 粒子模拟装置、粒子模拟方法以及粒子模拟程序
AU2017227323A AU2017227323B2 (en) 2016-03-03 2017-03-01 Particle simulation device, particle simulation method, and particle simulation program
GB1815732.1A GB2564326B (en) 2016-03-03 2017-03-01 Particle simulation device, particle simulation method, and particle simulation program
US16/080,875 US20200110911A1 (en) 2016-03-03 2017-03-01 Particle simulation device, particle simulation method, and particle simulation program

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
JP2016-041001 2016-03-03
JP2016041001A JP6671064B2 (ja) 2016-03-03 2016-03-03 粒子シミュレーション装置、粒子シミュレーション方法及び粒子シミュレーションプログラム

Publications (1)

Publication Number Publication Date
WO2017150626A1 true WO2017150626A1 (ja) 2017-09-08

Family

ID=59742999

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/JP2017/008168 WO2017150626A1 (ja) 2016-03-03 2017-03-01 粒子シミュレーション装置、粒子シミュレーション方法及び粒子シミュレーションプログラム

Country Status (6)

Country Link
US (1) US20200110911A1 (ja)
JP (1) JP6671064B2 (ja)
CN (1) CN108701275A (ja)
AU (1) AU2017227323B2 (ja)
GB (1) GB2564326B (ja)
WO (1) WO2017150626A1 (ja)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111104728A (zh) * 2019-10-15 2020-05-05 江汉大学 结构表面磨损的仿真预测方法和装置

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP7244757B2 (ja) * 2019-07-01 2023-03-23 富士通株式会社 情報処理装置、粒子シミュレーション方法および粒子シミュレーションシステム
JP2021110997A (ja) * 2020-01-07 2021-08-02 住友重機械工業株式会社 シミュレーション方法、シミュレーション装置、及びプログラム
CN111523272A (zh) * 2020-04-03 2020-08-11 北京三快在线科技有限公司 一种仿真测试的方法及装置
CN113705065A (zh) * 2021-08-20 2021-11-26 北京思朗科技有限责任公司 一种数据处理方法、装置、设备和可读存储介质
CN116187118A (zh) * 2022-11-25 2023-05-30 国家海洋环境预报中心 基于多gpu并行的多目标区域精细化海啸传播计算方法

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH11120165A (ja) * 1997-10-15 1999-04-30 Ricoh Co Ltd 粒子と壁との接触判定方法

Family Cites Families (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP4761865B2 (ja) * 2005-07-22 2011-08-31 株式会社日立製作所 空間シミュレーションプログラムおよび装置
JP5342761B2 (ja) * 2007-09-11 2013-11-13 プロメテック・ソフトウェア株式会社 粒子法に基づく流体シミュレーションの表面構築方法、そのプログラム、およびそのプログラムを格納した記憶媒体
JP5241468B2 (ja) * 2008-12-19 2013-07-17 住友重機械工業株式会社 シミュレーション方法及びプログラム
JP6261130B2 (ja) * 2014-06-04 2018-01-17 国立研究開発法人海洋研究開発機構 粒子シミュレーション装置、粒子シミュレーション方法及び粒子シミュレーションプログラム
CN104991999A (zh) * 2015-06-17 2015-10-21 大连理工大学 一种基于二维sph的溃坝洪水演进模拟方法

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH11120165A (ja) * 1997-10-15 1999-04-30 Ricoh Co Ltd 粒子と壁との接触判定方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
KOICHI RYU ET AL.: "Saikochiku Polygon-heki o Mochiita SPH no Kabe Kyokai Keisan Shuho ni Okeru Tansaku Han'i ni Tsuite", IPSJ SIG NOTES GRAPHICS TO CAD(CG) 2015-CG-160, vol. 2015 -CG, no. 5, 29 August 2015 (2015-08-29), pages 1 - 6 *
TAKASHI WATANABE ET AL.: "Study on the calculation procedure of wall boundary with large deformation by smoothed particle hydrodynamics method", JAPAN SOCIETY OF COMPUTATIONAL ENGINEERING AND SCIENCE RONBUNSHU, vol. 2013, 9 December 2013 (2013-12-09), XP055415484, ISSN: 1347-8826 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111104728A (zh) * 2019-10-15 2020-05-05 江汉大学 结构表面磨损的仿真预测方法和装置
CN111104728B (zh) * 2019-10-15 2023-04-25 江汉大学 结构表面磨损的仿真预测方法和装置

Also Published As

Publication number Publication date
US20200110911A1 (en) 2020-04-09
JP2017157076A (ja) 2017-09-07
AU2017227323A1 (en) 2018-09-13
JP6671064B2 (ja) 2020-03-25
CN108701275A (zh) 2018-10-23
AU2017227323B2 (en) 2019-10-03
GB2564326A (en) 2019-01-09
GB201815732D0 (en) 2018-11-14
GB2564326B (en) 2021-07-21

Similar Documents

Publication Publication Date Title
WO2017150626A1 (ja) 粒子シミュレーション装置、粒子シミュレーション方法及び粒子シミュレーションプログラム
JP6009075B2 (ja) 粒子流動のシミュレーションシステム及びその方法
Löhner Applied computational fluid dynamics techniques: an introduction based on finite element methods
JP5330300B2 (ja) 有限要素解析法における爆発シミュレーション
JP5800145B2 (ja) 解析装置、解析方法、解析プログラム及び記録媒体
EP2538359A2 (en) Boundary handling for particle-based simulation
Wachs Particle-scale computational approaches to model dry and saturated granular flows of non-Brownian, non-cohesive, and non-spherical rigid bodies
Valera et al. Modified algorithm for generating high volume fraction sphere packings
US20080162090A1 (en) System, methods, and computer readable media, for product design using t-spline deformation
US10503561B2 (en) Particle simulation apparatus and computer resource allocating method
Al Ibrahim et al. Particula: A simulator tool for computational rock physics of granular media
KR20200125459A (ko) 입자 거동 시뮬레이션 방법, 및 입자 거동 시뮬레이션 시스템
KR20090027142A (ko) 입자법에 기초한 유체 시뮬레이션의 표면 구축 방법, 그프로그램, 및 그 프로그램을 저장한 기억 매체
JP5113765B2 (ja) コンピュータシミュレーションおよび分析のための粒子への物体離散化
CN115618498A (zh) 一种飞行器跨流域流场的预测方法、装置、设备及介质
KR100588000B1 (ko) 유체 애니메이션에서의 자유경계 추적 장치 및 그 방법
JP5221603B2 (ja) 地盤変形解析装置、地盤変形解析方法、プログラム
Domínguez et al. Development of a new pre-processing tool for SPH models with complex geometries
Lee Developments in large scale discrete element simulations with polyhedral particles
JP5893333B2 (ja) 粒子挙動解析方法、粒子挙動解析装置、及び解析プログラム
JP5645120B2 (ja) 粒子状態計算装置及び粒子状態計算方法
Andersen et al. Post-processing in the material-point method
CN103377306B (zh) 在计算机辅助工程分析中提供颗粒连接模型的系统和方法
KR20110088043A (ko) 이동 그리드를 이용한 유체 시뮬레이션 시스템 및 방법
Playne et al. Benchmarking multi-GPU communication using the shallow water equations

Legal Events

Date Code Title Description
NENP Non-entry into the national phase

Ref country code: DE

ENP Entry into the national phase

Ref document number: 2017227323

Country of ref document: AU

Date of ref document: 20170301

Kind code of ref document: A

ENP Entry into the national phase

Ref document number: 201815732

Country of ref document: GB

Kind code of ref document: A

Free format text: PCT FILING DATE = 20170301

WWE Wipo information: entry into national phase

Ref document number: 1815732.1

Country of ref document: GB

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

Ref document number: 17760090

Country of ref document: EP

Kind code of ref document: A1

122 Ep: pct application non-entry in european phase

Ref document number: 17760090

Country of ref document: EP

Kind code of ref document: A1