US5848260A - Computer system for simulating physical processes - Google Patents
Computer system for simulating physical processes Download PDFInfo
- Publication number
- US5848260A US5848260A US08/566,032 US56603295A US5848260A US 5848260 A US5848260 A US 5848260A US 56603295 A US56603295 A US 56603295A US 5848260 A US5848260 A US 5848260A
- Authority
- US
- United States
- Prior art keywords
- elements
- voxels
- facet
- voxel
- sub
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Expired - Lifetime
Links
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T15/00—3D [Three Dimensional] image rendering
- G06T15/10—Geometric effects
- G06T15/40—Hidden part removal
-
- 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/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2111/00—Details relating to CAD techniques
- G06F2111/10—Numerical modelling
Definitions
- This invention relates to computer systems for simulating physical processes, e.g., fluid flow.
- the traditional lattice gas simulation assumes a limited number of particles at each lattice site, with the particles being represented by a short vector of bits.
- Each bit represents a particle moving in a particular direction.
- one bit in the vector might represent the presence (when set to 1) or absence (when set to 0) of a particle moving along a particular direction.
- Such a vector might have six bits, with, for example, the values 110000 indicating two particles moving in opposite directions along the X axis, and no particles moving along the Y and Z axes.
- a set of collision rules governs the behavior of collisions between particles at each site (e.g., a 110000 vector might become a 001100 vector, indicating that a collision between the two particles moving along the X axis produced two particles moving away along the Y axis).
- the rules are implemented by supplying the state vector to a lookup table, which performs a permutation on the bits (e.g., transforming the 110000 to 001100). Particles are then moved to adjoining sites (e.g., the two particles moving along the Y axis would be moved to neighboring sites to the left and right along the Y axis).
- Molvig et al. taught an improved lattice gas technique in which, among other things, many more bits were added to the state vector at each lattice site (e.g., 54 bits for subsonic flow) to provide variation in particle energy and movement direction, and collision rules involving subsets of the full state vector were employed.
- Molvig et al. taught an improved lattice gas technique in which, among other things, many more bits were added to the state vector at each lattice site (e.g., 54 bits for subsonic flow) to provide variation in particle energy and movement direction, and collision rules involving subsets of the full state vector were employed.
- the state vector instead of being a set of bits, was a set of integers (e.g., a set of eight bit bytes providing integers in the range of 0 to 255), each of which represented the number of particles in a given state.
- the system instead of being limited to a single particle moving in each direction at each momentum state, the system had the flexibility to model multiple particles moving in each direction at each momentum state.
- the Chen et al. system also provided a way of simulating the interaction between fluid particles and solid objects using a new "slip" technique that extended the simulation only to the outer surface of the boundary layer around a solid object, and not through it to the surface of the solid object.
- the collision rules governing interactions between particles and the surface allowed particles to retain tangential momentum. Instead of the particles being returned from the surface in the reverse tangential direction from which they arrived (as occurred with a prior "bounce back" technique), the particles were permitted to reflect from the surface and continue in the same tangential direction in which they arrived.
- the Chen et al. system specularly reflected particles with respect to the normal of the solid surface. A colliding particle left the surface at an opposite but approximately equal angle with respect to the normal.
- Chen et al.'s multi-particle model (in which more than one particle was allowed to have the same state at a point in the lattice) provided a simple and efficient way of modeling such slip collisions.
- the multi-particle model of Chen et al. accommodated arbitrary angular orientation of the solid boundary with respect to the lattice by allowing use of a weighted average of multiple outgoing states to assure that the average momentum of the outgoing particles was in a direction closely approximating true specular reflection. For example, from three to six different particle states (i.e., different directions and energies) could be used, each populated with enough outgoing particles that the net outgoing momentum was approximately equal to specular reflection of the group of incoming particles colliding with the surface.
- Chen et al. described techniques for preserving energy, mass and momentum normal to the solid boundary.
- Momentum normal to the solid boundary was preserved using a "pushing/pulling" technique that compared the overall incoming normal momentum to the overall outgoing normal momentum and recorded the normal surplus or deficit (i.e., the amount of normal momentum that had to be made up in some way that did not introduce artifacts into the simulation).
- Chen et al. then used a set of pushing/pulling rules to drive the normal surplus toward zero. Particles were moved from certain out-states to other out-states so that only normal momentum was affected.
- a cooling rule specified four states, two source states and two destination states.
- the two source states had the same momentum as the two destination states, but had a combined energy total that was higher than that of the two destination states. Moving a particle from both of the source states into both of the destination states decreased energy while conserving mass and momentum.
- Chen et al. used "dieting" rules to remove any surplus mass that accumulated as the result of one or more of the collision rules.
- a dieting rule specified two source states and one destination state. The two source states had a combined energy and momentum that was the same as the one destination state. Removing a particle from each of the source states and putting a particle into the one destination state caused the mass to be reduced while the energy and momentum were preserved.
- the invention features an improved technique for simulating the interaction between fluid particles and solid boundaries in a computer system for simulating physical processes.
- the surface of a boundary is represented as a collection of facets that are unrestricted by voxel boundaries, and may be arbitrarily oriented and sized relative to the lattice geometry of the voxels.
- the improved technique permits a continuously varying surface to be simulated with precision that is significantly greater than the precision of the lattice.
- a lattice gas is traditionally pictured as a set of discrete points distributed in a lattice structure that overlies a continuum. Particles of the lattice gas are point masses that propagate between lattice sites. The entire system is discretized in both space and time, and no physical meaning is assigned to non-integer space and time values. This model of a lattice gas has the virtue of being simple to envision and implement, but is not easily translated to a continual real physics system.
- An alternate physical picture of the lattice gas can be formed by assuming that the mass associated with each particle of the lattice gas is uniformly distributed throughout the volume of the voxel occupied by the particle. Propagation of particles may be viewed as movement of this distributed mass in the distance and direction indicated by the velocity of the particle. When particles are moving between adjacent fluid voxels, this distributed mass view of the lattice gas is equivalent to the traditional view of the lattice gas. However, the distributed mass view of the lattice gas offers the considerable advantage that it permits the particles of the lattice gas to interact with a surface that is arbitrarily oriented relative to the structure of the lattice.
- the invention features storing in a memory state vectors for multiple voxels and a representation for at least one surface that is sized and oriented independently of the size and orientation of the voxels.
- the state vectors include multiple entries that correspond to elements at particular momentum states of multiple possible momentum states at a voxel.
- Interaction operations that model interactions between elements of different momentum states are performed on the state vectors, and surface interaction operations that model interactions between a surface and elements at one or more voxels near the facet are performed on the representations of the surface.
- move operations that reflect movement of elements to new voxels are performed on the state vectors.
- the invention permits interactions with a surface to be simulated with increased resolution relative to the size of the voxels. This, in turn, permits highly accurate simulations of physical processes such as fluid flow.
- the entries of the state vectors may represents a density of elements per unit volume in a particular momentum state of a voxel. As discussed above, this distributive mass view of the state vectors offers considerable advantages.
- the entries of the state vectors may have integer values, and the interaction operations may be performed using integer values.
- the surface interaction operations may be performed using values representative of real number, such as floating point values.
- This hybrid approach permits the processing efficiency offered by integer operations to be obtained simultaneously with the structural flexibility offered by real numbers.
- exact conservation of mass, momentum and energy is enforced as simulated elements are moved from a surface facet to a voxel. Voxels that are intersected by a surface (and occupied by the region underlying the surface) may be viewed as partial voxels.
- a surface may be represented as a collection of facets that are sized and oriented independently of the size and orientation of the voxels and, in combination, represent at least one surface.
- the surface interaction operations model interactions between a facet and elements of at least one voxel near the facet.
- the elements may represent particles of a fluid and the facets may represent at least one surface over which the fluid flows.
- the surface interaction operations may include gathering elements from a first set of one or more voxels that interact with the facet, modelling interactions between the gathered elements and the facet to produce a set of surface interacted elements, and scattering the surface interacted elements to a second set of one or more voxels that interact with the facet.
- the surface interaction operations occur as a component of the operations by which elements are moved from voxel to voxel.
- the gather and scatter steps account for the arbitrary orientation of the facets. This permits the modelling step to be performed without regard to that orientation.
- the first set of voxels is identical to the second set of voxels.
- a factor indicative of the portion of the voxel that is not occupied by any facets (or material underlying a facet) is stored with the state vector for the facet.
- This factor permits operations to be performed on an occupied voxel without further regard to the degree to which the voxel is occupied.
- the factor for a voxel may be considered when gathering elements from the voxel.
- the elements obtained in the gather step for a particular momentum state may be identified as the elements occupying a region projected from the facet in the opposite direction of the velocity of the momentum state.
- the first set of voxels may include one or more voxels that overlap at least a portion of a region defined by a vector dot product of a vector representing a normal direction relative to the facet with a vector representing the velocity of the momentum state for which elements are being gathered.
- the first set of voxels may include all voxels that at least partially overlap a region projected from the facet and having a volume that equals the vector dot product multiplied by a surface area of the facet.
- the elements obtained in the gather step may also include elements scattered from one or more other facets.
- a surface includes non-convex features such as an interior corner, a region projected from a facet may be partially occupied by another facet. To account for this situation, elements are also gathered from appropriately oriented facets.
- interactions between the gathered elements and the facet are modelled so that the set of surface interacted elements satisfies one or more boundary conditions.
- the boundary conditions may include tangential momentum conditions. When the surface represented by a facet has no surface friction, the tangential momentum condition is satisfied when a combined momentum of the gathered elements that is tangential to the facet equals a combined momentum of the surface interacted elements that is tangential to the facet.
- the tangential momentum condition is satisfied when a difference between the combined tangential momentum of the gathered elements and the combined tangential momentum of the surface interacted elements is related to the level of surface friction of the surface (e.g., when the surface is completely non-slip, the combined tangential momentum of the surface interacted elements is zero).
- the boundary conditions may also include a normal momentum condition that is satisfied when a difference between a combined momentum of the gathered elements that is normal to the facet and a combined momentum of the surface interacted elements that is normal to the facet equals a pressure in a region defined by one or more voxels of the first set of voxels.
- the boundary conditions may also require conservation of mass and energy at the boundary. These conditions are satisfied when a combined mass of the gathered elements equals a combined mass of the surface interacted elements, and a combined energy of the gathered elements equals a combined energy of the surface interacted elements.
- a state vector includes integers representative of particular momentum states and the surface interaction operations are performed using values representative of real numbers, care must be taken to ensure that the integer values of the state vectors do not overflow or underflow when the surface interacted elements are scattered to the voxels. For this purpose, procedures that conserve the mass, momentum and energy of a voxel are employed whenever truncation is necessary to avoid overflow or underflow.
- the invention features simulating a physical process by storing in a memory state vectors for multiple voxels, with the state vectors including entries that correspond to particular momentum states of multiple possible momentum states at a voxel. Thereafter, operations are performed on the state vectors to model interactions between elements of different momentum states and movement of elements relative to the voxels, with at least some of the operations being performed using floating point numbers and some of the operations being performed using integers. Transitions between floating point and integer values may be made in a way that preserves physical properties of the physical process such as mass and energy.
- the entries of the state vectors may be integers representing a density of elements per unit volume in a particular momentum state.
- a representation of at least one surface that is sized and oriented independently of the size and orientation of the voxels may be stored in memory, and the operations performed using floating point numbers may be the operations that model interactions of the elements with the surface.
- FIG. 1 is a flow chart of a procedure followed by a physical process simulation system
- FIG. 2 is a perspective view of a microblock.
- FIGS. 3A and 3B are illustrations of lattice structures used by the system of FIG. 1.
- FIGS. 4 and 5 illustrate variable resolution techniques.
- FIG. 6 illustrates regions affected by a facet of a surface.
- FIG. 7 illustrates movement of particles from a voxel to a surface.
- FIG. 8 illustrates movement of particles from a surface to a surface.
- FIG. 9 is a flow chart of a procedure for performing surface dynamics.
- FIG. 10 illustrates an interface between voxels of different sizes.
- FIG. 11 is a flow chart of a procedure for simulating interactions with facets under variable resolution conditions.
- FIG. 12 is a block diagram of a functional unit of a physical process simulation system.
- FIG. 13 is a block diagram of a microdynamics unit of the system of FIG. 12.
- FIG. 14 is a block diagram of a single-voxel data path of the microdynamics unit of FIG. 13.
- a physical process simulation system operates according to a procedure 100 to simulate a physical process such as fluid flow.
- a simulation space is modelled as a collection of voxels (step 102).
- the simulation space is generated using a computer-aided-design (CAD) program.
- CAD computer-aided-design
- a CAD program could be used to draw an automobile positioned in a wind tunnel.
- data produced by the CAD program is processed to add a lattice structure having appropriate resolution and to account for objects and surfaces within the simulation space.
- Each lattice site, or voxel contains 54 states for subsonic mono-species simulations. The number of states is increased for transonic flows or multiple-species simulations.
- the resolution of the lattice may be selected based on the Reynolds number of the system being simulated.
- the Reynolds number is related to the viscosity (v) of the flow, the characteristic length (L) of an object in the flow, and the characteristic velocity (u) of the flow:
- the characteristic length of an object represents large scale features of the object. For example, if flow around an automobile were being simulated, the height of the automobile might be considered to be the characteristic length.
- the resolution of the simulation may be increased, or areas of increased resolution may be employed around the regions of interest.
- the dimensions of the voxels decrease as the resolution of the lattice increases.
- N i represents the number of elements, or particles, per unit volume in state i (i.e., the density of particles in state i) at a lattice site denoted by the three-dimensional vector x at a time t.
- N i represents the number of elements, or particles, per unit volume in state i (i.e., the density of particles in state i) at a lattice site denoted by the three-dimensional vector x at a time t.
- N i the number of particles
- the number of states is determined by the number of possible velocity vectors within each energy level.
- the velocity vectors consist of integer linear speeds in four-dimensional space: x, y, z and w.
- the fourth dimension, w is projected back onto three-dimensional space and thus does not indicate an actual velocity in the three-dimensional lattice.
- i ranges from 0 to 53.
- Each state represents a different velocity vector at a specific energy level (i.e., energy level zero, one or two).
- the velocity of each state is indicated with its "speed" in each of the four dimensions as follows:
- Energy level one states represent particles having a ⁇ 1 speed in two of the four dimensions and a zero speed in the other two dimensions
- Energy level two states represent particles having either a ⁇ 1 speed in all four dimensions, or a ⁇ 2 speed in one of the four dimensions and a zero speed in the other three dimensions.
- Each voxel (i.e., each lattice site) is represented by a state vector.
- the state vector completely defines the status of the voxel and includes 54 multi-bit entries, each of which corresponds to an integer value.
- the 54 entries correspond to the six rest states, 24 directional vectors at the first energy level and 24 directional vectors at the second energy level.
- Six rest states are employed to ensure that there are a sufficient number of rest "slots". Of course this same effect could be achieved by increasing the number of bits in the entry corresponding to the rest state in a 49 entry embodiment.
- the system offers substantially improved performance over systems that used single-bit entries to define voxel status.
- the system can produce Maxwell-Boltzmann statistics for the achieved equilibrium state vector.
- the voxels are grouped in 2 ⁇ 2 ⁇ 2 volumes that are called microblocks.
- the microblocks are organized to optimize parallel processing of the voxels as well as to minimize the overhead associated with the data structure.
- a short-hand notation for the voxels in the microblock is defined as N i (n), where n represents the relative position of the lattice site within the microblock and n ⁇ 0,1,2, . . . , 7 ⁇ .
- a microblock is illustrated in FIG. 2.
- a surface S (FIG. 3A) is represented in the simulation space (FIG. 3B) as a collection of facets F.sub. ⁇ :
- each facet F.sub. ⁇ has a unit normal (n.sub. ⁇ ), a surface area (A.sub. ⁇ ), a center location (x.sub. ⁇ ), and facet distribution function (N i ( ⁇ )) that describes the surface dynamic properties of the facet.
- different levels of resolution may be used in different regions of the simulation space to improve processing efficiency.
- the region 150 around an object 155 is of the most interest and is therefore simulated with the highest resolution.
- decreasing levels of resolution i.e., expanded voxel volumes
- regions 160, 165 that are spaced at increasing distances from the object 155.
- a lower level of resolution may be used to simulate a region 170 around less significant features of an object 175 while the highest level of resolution is used to simulate regions 180 around the most significant features (e.g., the leading and trailing surfaces) of the object 175.
- Outlying regions 185 are simulated using the lowest level of resolution and the largest voxels. Techniques for processing interactions between differently sized voxels are discussed by Molvig et al. in U.S. Pat. No. 5,377,129, which is incorporated herein by reference, at col. 18, line 58 to col. 28, line 21. Techniques for processing interactions between differently sized voxels and the facets of a surface are discussed below.
- Voxels affected by one or more facets are identified (step 104).
- Voxels may be affected by facets in a number of ways.
- voxels that are intersected by facets are affected in that each of these voxels has a reduced volume relative to non-intersected voxels. This occurs because a facet, and material underlying the surface represented by the facet, occupies a portion of the voxel.
- a fractional factor P .sub. ⁇ (x) indicates the portion of the voxel that is unaffected by the facet (i.e., the portion that can be occupied by a fluid or other materials for which flow is being simulated). For non-intersected voxels, P.sub. ⁇ (x) equals one.
- Voxels that interact with one or more facets by transferring particles to the facet or receiving particles from the facet are also identified as voxels affected by the facets. All voxels that are intersected by a facet will include at least one state that receives particles from the facet and at least one state that transfers particles to the facet. In most cases, additional voxels also will include such states.
- a facet F.sub. ⁇ receives particles from, or transfers particles to, a region defined by a parallelepiped G i ⁇ having a height defined by the magnitude of the vector dot product of the velocity vector c i and the unit normal n.sub. ⁇ of the facet (
- the facet F.sub. ⁇ receives particles from the volume V i ⁇ when the velocity vector of the state is directed toward the facet (c i ⁇ n.sub. ⁇ ⁇ 0), and transfers particles to the region when the velocity vector of the state is directed away from the facet (c i ⁇ n.sub. ⁇ >0).
- this expression must be modified when another facet occupies a portion of the parallelepiped G i ⁇ , a condition that could occur in the vicinity of non-convex features such as interior corners.
- the parallelepiped G i ⁇ of a facet F.sub. ⁇ may overlap portions or all of multiple voxels.
- the number of voxels or portions thereof is dependent on the size of the facet relative to the size of the voxels, the energy of the state, and the orientation of the facet relative to the lattice structure.
- the number of affected voxels increases with the size of the facet. Accordingly, the size of the facet, as noted above, is typically selected to be on the order of or smaller than the size of the voxels located near the facet.
- V i ⁇ (x) The portion of a voxel N(x) overlapped by a parallelepiped G i ⁇ is defined as V i ⁇ (x).
- the number of state i particles that move between a voxel N(x) and a facet F.sub. ⁇ equals the density of state i particles in the voxel multiplied by the volume of the region of overlap with the voxel:
- a timer is initialized to begin the simulation (step 106).
- movement of particles from voxel to voxel is simulated by an advection stage (steps 108-116) that accounts for interactions of the particles with surface facets.
- a collision stage (step 118) simulates the interaction of particles within each voxel.
- the timer is incremented (step 120). If the incremented timer does not indicate that the simulation is complete (step 122), the advection and collision stages (steps 108-120) are repeated. If the incremented timer indicates that the simulation is complete (step 122), results of the simulation are stored and/or displayed (step 124).
- each facet must meet four boundary conditions.
- the combined mass of particles received by a facet must equal the combined mass of particles transferred by the facet (i.e., the net mass flux to the facet must equal zero).
- the combined energy of particles received by a facet must equal the combined energy of particles transferred by the facet (i.e., the net energy flux to the facet must equal zero).
- the other two boundary conditions are related to the net momentum of particles interacting with a facet.
- a slip surface For a surface with no skin friction, referred to herein as a slip surface, the net tangential momentum flux must equal zero and the net normal momentum flux must equal the local pressure at the facet.
- the components of the combined received and transferred momentums that are perpendicular to the normal n.sub. ⁇ of the facet i.e., the tangential components
- the difference between the components of the combined received and transferred momentums that are parallel to the normal n.sub. ⁇ of the facet i.e., the normal components
- friction of the surface reduces the combined tangential momentum of particles transferred by the facet relative to the combined tangential momentum of particles received by the facet by a factor that is related to the amount of friction.
- step 108 particles are gathered from the voxels and provided to the facets (step 108).
- the number of state i particles that move between a voxel N(x) and a facet F.sub. ⁇ is:
- particles are moved between facets (step 110). If the parallelepiped G i ⁇ for an incoming state (c i ⁇ n.sub. ⁇ ⁇ 0) of a facet F.sub. ⁇ is intersected by another facet F.sub. ⁇ , then a portion of the state i particles received by the facet F.sub. ⁇ will come from the F.sub. ⁇ . In particular, facet F.sub. ⁇ will receive a portion of the state i particles produced by facet F.sub. ⁇ during the previous time increment. This relationship is illustrated in FIG.
- ⁇ i ( ⁇ ) is a measure of the state i particles produced by the facet F.sub. ⁇ during the previous time increment. From this, for each state i directed toward a facet F.sub. ⁇ (c i ⁇ n.sub. ⁇ 0), the number of particles provided to the facet F.sub. ⁇ by the other facets is: ##EQU4## and the total flux of state i particles into the facet is: ##EQU5##
- the state vector V( ⁇ ) for the facet also referred to as a facet distribution function, has 54 entries corresponding to the 54 entries of the voxel state vectors.
- the input states of the facet distribution function N( ⁇ ) are set equal to the flux of particles into those states divided by the volume V i ⁇ :
- the facet distribution function is a simulation tool for generating the output flux from a facet, and is not necessarily representative of actual particles. To generate an accurate output flux, values are assigned to the other states of the distribution function. Outward states are populated using the techniques described above for populating the inward states:
- ⁇ iOTHER ( ⁇ ) is determined using the technique described above for generating ⁇ iIN ( ⁇ ).
- ⁇ iOTHER ( ⁇ ) is generating using values of ⁇ iOUT ( ⁇ ) from the previous time step so that:
- N i ( ⁇ ) for parallel states is determined as the limit of N i ( ⁇ ) as V i ⁇ and V i ⁇ (x) approach zero.
- states having zero velocity i.e., rest states and states (0, 0, 0, 2) and (0, 0, 0,-2) are initialized at the beginning of the simulation based on initial conditions for temperature and pressure. These values are then adjusted over time.
- step 112 surface dynamics are performed for each facet to satisfy the four boundary conditions discussed above (step 112).
- a procedure for performing surface dynamics for a facet is illustrated in FIG. 9. Initially, the combined momentum normal to the facet F x is determined (step 905) by determining the combined momentum P( ⁇ ) of the particles at the facet as: ##EQU6## From this, the normal momentum P n ( ⁇ ) is determined as:
- An outgoing flux distribution for the facet F.sub. ⁇ is then determined (step 920) based on the incoming flux distribution and the Boltzmann distribution.
- the difference between the incoming flux distribution ⁇ i ( ⁇ ) and the Boltzmann distribution is determined as:
- i* is the state having a direction opposite to state i. For example, if state i is (1, 1, 0, 0) , then state i* is (-1, -1, 0, 0).
- the outgoing flux distribution may be further refined to: ##EQU7## for n.sub. ⁇ ⁇ c i >0 and where K.sub. ⁇ is a function of skin friction, t 1 ⁇ is a first tangential vector that is perpendicular to n.sub. ⁇ , t 2a is a second tangential vector that is perpendicular to both n.sub. ⁇ and t 1 ⁇ , and ⁇ N j ,1 and ⁇ N j ,2 are distribution functions corresponding to the energy (j) of the state i and the indicated tangential vector.
- the distribution functions are determined according to: ##EQU8## where j equals 1 for energy level 1 states and 2 for energy level 2 states.
- each line of the equation for ⁇ iOUT ( ⁇ ) are as follows.
- the terms on the first line enforce the normal momentum flux boundary condition to the extent that collisions have been effective in producing a Boltzmann distribution, but include a tangential momentum flux anomaly.
- the terms on the third and fourth lines correct for this anomaly, which may arise due to discreteness effects or non-Boltzmann structure due to insufficient collisions.
- the term on the second line adds a specified amount of skin friction to enforce a desired change in tangential momentum flux on the surface. Note that all terms involving vector manipulations are geometric factors that may be calculated prior to beginning the simulation.
- ⁇ is the density of the facet distribution: ##EQU9##
- T of the facet distribution is then determined as: ##EQU10## where j equals 0, 1 or 2 and corresponds to the energy of the state i.
- d 0 is the number of rest states (i.e., six).
- the difference between the input energy and the output energy is measured for each energy level j as: ##EQU14## where the index j denotes the energy of the state i. This energy difference is then used to generate a difference term: ##EQU15## for c ji ⁇ n.sub. ⁇ >0. This difference term is used to modify the outgoing flux so that the flux becomes:
- particles are moved between voxels along the three-dimensional rectilinear lattice (step 114).
- This voxel to voxel movement is the only movement operation performed on voxels that do not interact with the facets (i.e., voxels that are not located near a surface).
- voxels that are not located near enough to a surface to interact with the surface constitute a large majority of the voxels.
- Each of the separate states represents particles moving along the lattice with integer speeds in each of the four dimensions: x, y, z and w.
- the integer speeds include: 0, ⁇ 1, and ⁇ 2.
- the sign of the speed indicates the direction in which a particle is moving along the corresponding axis.
- the move operation is computationally quite simple.
- the entire population of a state is moved from its current voxel to its destination voxel during every time increment.
- the particles of the destination voxel are moved from that voxel to their own destination voxels.
- an energy level 1 particle that is moving in the +1x and +1y direction (1, 1, 0, 0) is moved from its current voxel to one that is +1 over in the x direction and +1 up in the y direction.
- the particle ends up at its destination voxel with the same state it had before the move (1,1,0,0). Interactions within the voxel will likely change the particle count for that state based on local interactions with other particles and surfaces. If not, the particle will continue to move along the lattice at the same speed and direction.
- particles in the stopped states do not get moved.
- the particles only get moved in three dimensions.
- Non-zero values in the w dimension of a state do not affect the determination of a lattice site to which the state's particles are to be moved. For example, an energy level one particle with a -1z and a +1w speed (0,0,-1,1) and an energy level one particle with a -1z and a -1w (0,0,-1,-1) speed would both move to a site that is -1 away in the z dimension.
- V i ⁇ (x) and P.sub. ⁇ (x) may have non-integer values. This can result in one or more fractional particles being transferred to a facet. Transfer of such fractional particles to a facet results in fractional particles remaining in the voxels. These fractional particles are transferred to a voxel occupied by the facet. For example, referring to FIG.
- step 108 when a portion 700 of the state i particles for a voxel 705 is moved to a facet 710 (step 108), the remaining portion 715 is moved to a voxel 720 in which the facet 710 is located and from which particles of state i are directed to the facet 710.
- the state population equalled 25 and V i ⁇ (x) equalled 0.25 (i.e., a quarter of the voxel intersects the parallelepiped G i ⁇ )
- 6.25 particles would be moved to the facet F.sub. ⁇ and 18.75 particles would be moved to the voxel occupied by the facet F.sub. ⁇ .
- N(f) is the source voxel. Because N i (f) may be a real number, N i (f) is temporarily stored in a buffer as a floating point number or a scaled integer until the advection stage is completed.
- N i (f) is determined as: ##EQU17## where scale is a constant having a value of 2 16 and N si (x) and V si ⁇ (x) are, respectively, N i (x) and V i ⁇ (x) multiplied by scale.
- the random number rand takes on values between zero and one prior to scaling, and is then scaled up by scale to values between 0 and 2 16 -1.
- the operator floor then produces an integer value. In particular, floor returns the largest integer that is less than or equal to its argument. For example, floor(2.3) equals 2 and floor(-2.3) equals -3.
- step 116 the outgoing particles from each facet are scattered to the voxels.
- this step is the reverse of the gather step by which particles were moved from the voxels to the facets.
- the number of state i particles that move from a facet F.sub. ⁇ to a voxel N(x) is: ##EQU18## where P.sub. ⁇ (x) accounts for the volume reduction of partial voxels.
- t s- is the time just before the scatter step
- t s + is the time just after the scatter step
- N i (x, t s- ) is a floating point number.
- the random number rand takes on values between zero and one prior to scaling, and is then scaled up by scale to values between 0 and 2 16 -1. Because the values are scaled up by scale, the division by scale is merely a matter of shifting the result of the multiplication and addition operations.
- an amount of mass equal to the value gained (due to underflow) or lost (due to overflow) is added back to randomly (or sequentially) selected states having the same energy and that are not themselves subject to overflow or underflow.
- the additional momentum resulting from this addition of mass and energy is accumulated and added to the momentum from the truncation.
- both mass and energy are corrected when the mass counter reaches zero.
- the momentum is corrected using pushing/pulling techniques until the momentum accumulator is returned to zero.
- fluid dynamics are performed (step 118). This step may be referred to as microdynamics or intravoxel operations. Similarly, the advection procedure may be referred to as intervoxel operations.
- the microdynamics operations described below may also be used to collide particles at a facet to produce a Boltzmann distribution.
- the microdynamics operations simulate the set of physical interactions that occur within a voxel by sequentially applying a series of binary, trinary, or "n"ary interaction rules (also referred to as collision rules) to the state vector of a voxel during each time increment.
- n binary, trinary, or "n”ary interaction rules
- collision rules there are 276 collision rules, each of which is implemented in a separate collision stage. Since the collision rules are applied sequentially, the result of their application has a complex nonlinear dependence on the initial state values that cannot be written simply in terms of those values.
- the collision rules can be implemented as binary, trinary, or "n”ary collision events, for simplicity and ease of understanding the following discussion will refer primarily to binary collision events.
- Non-energy exchanging collisions allow for particles to collide with each other, thus changing their velocity. Because the state occupied by a particle determines the particle's velocity, a change in a particle's velocity is accomplished by moving that particle into a different state.
- a binary collision rule describes a self collision with reference to two input states and two output states, and simulates the collision of one or more particles from each of the input states by moving the particles to the output states. Because the collision rules are bi-directional, the input and output states are determined at the time the collision takes place depending on the states' populations.
- Energy exchanging collisions differ from self collisions in that the two outgoing particles are at different energy levels than the two incoming particles. For subsonic flows there are only three energy levels: zero (stopped), one and two. To conserve energy, the only possible energy exchanging collisions occur when one pair includes two energy level one particles and the other pair includes an energy level two particle and a stopped particle. Energy exchanging collisions do not happen at the same rate in both directions. Rather, they happen at a forward rate from the energy level one states and a backward rate from the energy two and stopped states. As will be discussed in more detail below, these collision rates are dependent on the temperature of the system.
- Each collision rule whether directed to self or energy exchanging collisions, operates on a subset of the states of the state vector for a single voxel or facet and simulates the net effect of collisions between particles in those states.
- a binary collision rule operates on four states (i, j, k, l) and simulates the net effect of collisions between particles in the first two states (i, j) and collisions between particles in the second two states (k, l).
- a binary collision rule can be expressed as: ##EQU20## where ⁇ denotes the physically infinitesimal time required to perform the collision rule.
- a trinary collision rule can be expressed as: ##EQU21##
- a state change resulting from application of a collision rule actually reflects a net state change.
- this movement might actually represent collisions that resulted in, for example, fifty six particles from each of states i and j moving to each of states k and l, and fifty two particles from each of states k and l moving to each of states i and j.
- each collision rule must conserve mass, momentum and energy. Because the rules simply move particles from one state to another, and do not create or destroy particles, they necessarily conserve mass. To ensure that a rule conserves momentum, the states affected by the rule are selected so that:
- the states affected by the rule are selected so that:
- the states i, j, k, and l are selected for each rule so that a pair of particles from states i and j has the same total momentum and energy as a pair of particles from states k and l. Because the incoming and outgoing pairs of states must always conserve mass, momentum and energy, not all possible quartets within the 54 states correspond to a "legal" collision rule.
- N scatt is determined by applying a collision operator to the states affected by the rule.
- the arithmetic/computational form of the collision operator determines the practicality of implementing the operator in a simulation system.
- the statistical properties of the collision operator determine the extent to which accurate behavior can be achieved during a simulation.
- the preferred collision operator is one for which N scatt approaches zero as the number of particles in each state affected by a rule approaches its Boltzmann equilibrium value.
- the preferred operator is the "multilinear rule", which has the following form:
- a j is a collision coefficient and R.sub. ⁇ and R b are the forward and backward collision rates for energy exchanging collisions.
- the coefficients A j depend only on local temperature and density, while R.sub. ⁇ and R b depend only on local temperature. Since temperature and density do not change due to collisions, A j , R.sub. ⁇ and R b are constants for a particular time increment, and can be determined prior to all collisions.
- states i and j are the energy level one states and states k and l are the energy level two and rest states so that R.sub. ⁇ /R b represents the rate of collisions from energy level one particles to energy level two and rest particles.
- R.sub. ⁇ /R b is defined as: ##EQU23## where d 0 is the number of stopped particle states and equals six in the described embodiment.
- the temperature of the fluid is not necessarily constant over the length of a simulation, especially for simulations involving heat transfer, and A j and R.sub. ⁇ /R b must be updated dynamically during the simulation to reflect changes in the local temperature.
- the temperature range supported for subsonic flows is between 1/3 and 2/3.
- R.sub. ⁇ /R b equals one.
- i, j, k and l are selected so that the combined momentum of states i and j and of states k and l are two in the x dimension and zero in the y, z and w dimensions, and so that each state is an energy level one state.
- ⁇ is 1656 and the temperature is 1/2
- a l using the expression for ⁇ N> l .sup.(eq) so that ⁇ N> l .sup.(eq) equals ⁇ /48 and A l , equals 12/ ⁇ or 1/138.
- the collision rule were applied using floating point numbers the resulting state populations would have the following values:
- N scatt is maintained as an integer value.
- N scatt is determined as:
- N scatt is determined as:
- R.sub. ⁇ ' and R b ' are, respectively, R.sub. ⁇ and R b multiplied by A.sub. ⁇ s, which is A.sub. ⁇ multiplied by scale.
- Rand is employed to prevent the truncation that results from the floor operation from introducing statistical bias into the system. Use of rand ensures that this operation, which forces N scatt to take on integer values, will not statistically bias N scatt in a particular direction.
- the values of N i , N j , N k and N l for the example provided above would be set to (25, 40, 53, 20) with 56.522% probability and (26, 41, 52, 19) with 43.478% probability.
- the mean values for N i , N j , N k and N l would be (25.43478, 40.43478, 52.56522, 19.56522).
- the multilinear collision operator drives the simulated system to Boltzmann equilibrium. This means that repeated application of the collision rules drives the system to an equilibrium value where, for a given set of macroscopic conditions, the individual state populations, N i , take on known values as determined by the Boltzmann distribution: ##EQU25## where the factors in the exponent represent the conserved invariants of mass (m), momentum (mc i ) and energy (mc i 2 ), it is assumed that no additional (“spurious") invariants occur, and r i represents a weighting factor generated from R.sub. ⁇ and R b . Attainment of this distribution after each application of the collision rules will result in hydrodynamic behavior to some accuracy.
- N i N i must be a sum of summational invariants of the collision process (i.e. quantities whose sum over particles are conserved) in the following form: ##EQU26##
- I i .sup. ⁇ is the ⁇ th invariant associated with state i, (e.g., a component of the momentum, m i c i , in the equation for N i above)
- A.sub. ⁇ are the coefficients associated with their corresponding invariants (e.g., ⁇ , ⁇ and ⁇ in the equation for N i above).
- variable resolution employs voxels of different sizes, hereinafter referred to as coarse voxels 1000 and fine voxels 1005.
- coarse voxels 1000 and fine voxels 1005.
- fine voxels 1005. The interface between regions of coarse and fine voxels is referred to as a variable resolution (VR) interface 1010.
- facets may interact with voxels on both sides of the VR interface. These facets are classified as VR interface facets 1015 (F.sub. ⁇ IC) or VR fine facets 1020 (F.sub. ⁇ IF).
- a VR interface facet 1015 is a facet positioned on the coarse side of the VR interface and having a coarse parallelepiped 1025 extending into a fine voxel.
- a coarse parallelepiped is one for which c i is dimensioned according to the dimensions of a coarse voxel
- a fine parallelepiped is one for which c i is dimensioned according to the dimensions of a fine voxel.
- a VR fine facet 1020 is a facet positioned on the fine side of the VR interface and having a fine parallelepiped 1030 extending into a coarse voxel. Processing related to interface facets may also involve interactions with coarse facets 1035 (F.sub. ⁇ C) and fine facets 1040 (F.sub. ⁇ F).
- VR facets For both types of VR facets, surface dynamics are performed at the fine scale, and operate as described above. However, VR facets differ from other facets with respect to the way in which particles advect to and from the VR facets.
- Interactions with VR facets are handled using a variable resolution procedure 1100 illustrated in FIG. 11. Most steps of this procedure are carried out using the comparable steps discussed above for interactions with non-VR facets.
- the procedure 1100 is performed during a coarse time step (i.e., a time period corresponding to a coarse voxel) that includes two phases that each correspond to a fine time step.
- the facet surface dynamics are performed during each fine time step.
- a VR interface facet F.sub. ⁇ IC is considered as two identically sized and oriented fine facets that are referred to, respectively, as a black facet F.sub. ⁇ ICb and a red facet F.sub. ⁇ ICr.
- the black facet F.sub. ⁇ ICb is associated with the first fine time step within a coarse time step while the red facet F.sub. ⁇ ICr is associated with the second fine time step within a coarse time step.
- particles are moved (advected) between facets by a first surface-to-surface advection stage (step 1102).
- Particles are moved from black facets F.sub. ⁇ ICb to coarse facets F.sub. ⁇ C with a weighting factor of V.sub. ⁇ that corresponds to the volume of the unblocked portion of the coarse parallelepiped (FIG. 10, 1025) that extends from a facet F.sub. ⁇ and that lies behind a facet F.sub. ⁇ less the unblocked portion of the fine parallelepiped (FIG. 10, 1045) that extends from the facet F.sub. ⁇ and that lies behind the facet F.sub. ⁇ .
- the magnitude of c i for a fine voxel is one half the magnitude of c i for a coarse voxel.
- the volume of a parallelepiped for a facet F.sub. ⁇ is defined as:
- the volume of a fine parallelepiped corresponding to a facet is one half the volume of the corresponding coarse parallelepiped for the facet.
- Particles are moved from coarse facets F.sub. ⁇ C to black facets F.sub. ⁇ ICb with a weighting factor of V.sub. ⁇ that corresponds to the volume of the unblocked portion of the fine parallelepiped that extends from a facet F.sub. ⁇ and that lies behind a facet F.sub. ⁇ .
- Particles are moved from red facets F.sub. ⁇ ICr to coarse facets F.sub. ⁇ C with a weighting factor of V.sub. ⁇ , and from coarse facets F.sub. ⁇ C to red facets F.sub. ⁇ ICr with a weighting factor of V.sub. ⁇ .
- Particles are moved from red facets F.sub. ⁇ ICr to black facets F.sub. ⁇ ICb with a weighting factor V.sub. ⁇ .
- black-to-red advections do not occur.
- black and red facets represent consecutive time steps, black-to-black advections (or red-to-red advections) never occur.
- particles in this stage are moved from red facets F.sub. ⁇ ICr to fine facets F.sub. ⁇ IF or F.sub. ⁇ F with a weighting factor V.sub. ⁇ , and from fine facets F.sub. ⁇ IF or F.sub. ⁇ F to black facets F.sub. ⁇ ICb with the same weighting factor.
- particles are moved from fine facets F.sub. ⁇ IF or F.sub. ⁇ F other fine facets F.sub. ⁇ IF or F.sub. ⁇ F with the same weighting factor and from coarse facets F.sub. ⁇ C to other coarse facets F C with a weighting factor V C ⁇ that corresponds to the volume of the unblocked portion of the coarse parallelepiped that extends from a facet F.sub. ⁇ and that lies behind a facet F.sub. ⁇ .
- particles are gathered from the voxels in a first gather stage (steps 1104-1110). Particles are gathered for fine facets F.sub. ⁇ F from fine voxels using fine parallelepipeds (step 1104), and for coarse facets F.sub. ⁇ C from coarse voxels using coarse parallelepipeds (step 1106). Particles are then gathered for black facets F.sub. ⁇ IRb and for VR fine facets F.sub. ⁇ IF from both coarse and fine voxels using fine parallelepipeds (step 1108). Finally, particles are gathered for red facets F.sub. ⁇ IRr from coarse voxels using the differences between coarse parallelepipeds and fine paralllelepipeds (step 1110).
- coarse voxels that interact with fine voxels or VR facets are exploded into a collection of fine voxels (step 1112).
- the states of a coarse voxel that will transmit particles to a fine voxel within a single coarse time step are exploded.
- the appropriate states of a coarse voxel that is not intersected by a facet are exploded into eight fine voxels oriented like the microblock of FIG. 2.
- the appropriate states of coarse voxel that is intersected by one or more facets are exploded into a collection of complete and/or partial fine voxels corresponding to the portion of the coarse voxel that is not intersected by any facets.
- the particle densities N i (x) for a coarse voxel and the fine voxels resulting from the explosion thereof are equal, but the fine voxels may have fractional factors P.sub. ⁇ that differ from the fractional factor of the coarse voxel and from the fractional factors of the other fine voxels.
- particles are moved between fine voxels (step 1118) including actual fine voxels and fine voxels resulting from the explosion of coarse voxels. Once the particles have been moved, particles are scattered from the fine facets F.sub. ⁇ IF and F.sub. ⁇ F to the fine voxels (step 1120).
- Particles are also scattered from the black facets F.sub. ⁇ ICb to the fine voxels (including the fine voxels that result from exploding a coarse voxel) (step 1122). Particles are scattered to a fine voxel if the voxel would have received particles at that time absent the presence of a surface.
- particles are scattered to a voxel N(x) when the voxel is an actual fine voxel (as opposed to a fine voxel resulting from the explosion of a coarse voxel), when a voxel N(x+c i ) that is one velocity unit beyond the voxel N(x) is an actual fine voxel, or when the voxel N(x+c i ) that is one velocity unit beyond the voxel N(x) is a fine voxel resulting from the explosion of a coarse voxel.
- the first fine time step is completed by performing fluid dynamics on the fine voxels (step 1124).
- the voxels for which fluid dynamics are performed do not include the fine voxels that result from exploding a coarse voxel (step 1112).
- the procedure 1100 implements similar steps during the second fine time step. Initially, particles are moved between surfaces in a second surface-to-surface advection stage (step 1126). Particles are advected from black facets to red facets, from black facets to fine facets, from fine facets to red facets, and from fine facets to fine facets.
- particles are advected between surfaces particles are gathered from the voxels in a second gather stage (steps 1128-1130).
- Particles are gathered for red facets F.sub. ⁇ IRr from fine voxels using fine parallelepipeds (step 1128).
- Particles also are gathered for fine facets F.sub. ⁇ F and F.sub. ⁇ IF from fine voxels using fine parallelepipeds (step 1130).
- step 1132 surface dynamics are performed for the fine facets F.sub. ⁇ IF and F.sub. ⁇ F (step 1132), for the coarse facets F.sub. ⁇ C (step 1134), and for the red facets F.sub. ⁇ ICr (step 1136) as discussed above.
- particles are moved between voxels using fine resolution (step 1138) so that particles are moved to and from fine voxels and fine voxels representative of coarse voxels.
- Particles are then moved between voxels using coarse resolution (step 1140) so that particles are moved to and from coarse voxels.
- particles are scattered from the facets to the voxels while the fine voxels that represent coarse voxels (i.e., the fine voxels resulting from exploding coarse voxels) are coalesced into coarse voxels (step 1142).
- particles are scattered from coarse facets to coarse voxels using coarse parallelepipeds, from fine facets to fine voxels using fine parallelepipeds, from red facets to fine or coarse voxels using fine parallelepipeds, and from black facets to coarse voxels using the differences between coarse parallelepipeds and find parallelepipeds.
- fluid dynamics are performed for the fine voxels and the coarse voxels (step 1144).
- ASIC 12 application specific integrated circuit
- memory bank 14 the system described above can be implemented using a functional unit 10 that includes an application specific integrated circuit ("ASIC") 12 coupled to a memory bank 14 and a rule store 16.
- ASIC 12 is designed so that hundreds, or even thousands, of functional units 10 can be connected together and operated in parallel for improved performance.
- ASIC 12 includes a processor control unit (“PCU”) 18 that itself includes a reduced instruction set computer (“RISC”) processor 20, an instruction cache 22, and RISC interface logic 24.
- PCU processor control unit
- RISC reduced instruction set computer
- Processor control unit 18 acts as a central controller and issues vector-like commands to the various functional units.
- the functional units act as co-processors, taking an instruction from the processor control unit 18 and executing their own set of predefined sequences of instructions to complete the task.
- the instructions for processor control unit 18 are downloaded into memory bank 14 and are executed from there by processor control unit 18. These instructions control processor control unit 18 and cause it to issue commands to the other functional units in ASIC 12.
- memory bank 14 stores particle states for each voxel assigned to the ASIC 12 along with the associated data structure and overhead information.
- the amount of storage provided by memory bank 14 dictates the amount of resolution that ASIC 12 can be assigned as its part of the total simulation volume.
- Memory bank 14 is constructed from a set of standard commodity synchronous DRAMs, which have been selected over standard asynchronous DRAMS due to their pipelined access mode and high data bandwidth capabilities.
- a microdynamics unit which includes microdynamics control logic 26 and a microdynamics data path 28, is the functional unit that is responsible for processing the particle dynamics. This includes particle interactions with each other as well as with surfaces.
- Microdynamics data path 28 includes the hardware necessary to perform the computations on the voxel data.
- Microdynamics control logic 26 sets up the data path to perform the correct operation and selects the proper data to be processed.
- Microdynamics control logic 26 sets up the data path and selects the data based on rule sets retrieved from rule store 16.
- Rule store 16 is constructed out of SRAMs to provide a fast access time and the ability to load the memory as the microdynamics rules are changed.
- a multi-port RAM 30 acts as a register file for the particle states and the associated data that the microdynamics data path 28 is currently processing.
- RAM 30, which is a six-port static RAM, has 256 64-bit entries and supplies two read ports and two write ports to microdynamics data path 28.
- RAM 30 provides a read port and a write port to a microblock unit 32. This third set of ports is used for concurrent data movement between RAM 30 and a memory access unit 34 as well as for statistics.
- Memory access unit 34 is the central source and destination of particle data transfers within ASIC 12, and controls memory bank 14.
- microdynamics data path 28 processes eight voxels in parallel, by providing eight identical copies of the eight-bit data path 29 necessary to process one voxel.
- Microdynamics control logic 26 is constructed as a SIMD (single instruction multiple data) controller, where all eight sites are set up to execute the same operation on their independent voxel data.
- a bus interface unit (“BIU") 36 connects ASIC 12 to external processors such as, for example, a general purpose "host” processor used for system control. (The host processor, or another external processor, would be used to implement the floating point operations employed when modelling the interactions of particles with faceted surfaces.)
- the main functions of BIU 36 are to provide a path for the system's host processor to provide initialization information to ASIC 12, and to send statistical information from ASIC 12 to the host processor.
- BIU 36 provides a path to memory access unit 34 so that memory bank 14 is accessible to the host processor.
- An advection unit (“AU") 38 provides the data path and control to enables particles to advect (move) along the cubic lattice, including support for those particles that are moving to a destination voxel that is owned by another ASIC 12. To transmit particles to a voxel owned by another ASIC 12, AU 38 includes communications ports.
- Microblock unit 32 is responsible for calculating statistics on the data being processed in the MDU. After the voxel data has been processed in the MDU, microblock unit 32 uses the voxel data to generate the statistics. The data is then sent back to memory bank 14 through memory access unit 34.
- the RISC processor 20 is capable of doing the standard set of 32-bit integer arithmetic and logical operations. In effect, the RISC processor can be coded to execute a program in memory as any standard processor could.
- a memory map has been established for ASIC 12 that essentially creates three segments in the processors 4 GB address space:
- the logic in the PCU decodes the address issued by the RISC processor and determines where it should be sent. Local memory requests are sent to the MAU. Memory mapped register reads or writes which includes functional unit instructions are sent to the appropriate functional unit. Host memory references are sent to the BIU which will forward the requests to the host system.
- the MAU arbitrates between the requests from these four units and allows single and multiple word transfers.
- the MAU maintains a DMA engine to support the burst transfers.
- the unit Upon receiving a grant from the MAU, the unit connects to the MAU and drives or receives data until its transfer is complete.
- the RISC processor initiates the transfer of most of the data packets between the various functional units and the MAU.
- the function units, or co-processors in this model have no knowledge of system or ASIC level transactions.
- the units are given data and tasks to perform by the RISC engine and they process that data locally.
- the processor instructs the microblock unit to fetch a specific block of data from the memory controlled by the MAU and place it into the multi-port RAM.
- the RISC engine then issues an instruction to the microdynamics unit to process that data.
- the RISC engine instructs the microblock unit to retrieve the block of data from the multi-port RAM and put it back in a specific location in memory.
- the co-processor units are capable of processing routines that take from tens to hundreds to even thousands of cycles based on a single command from the RISC engine and a block of data held in local storage. This processing length allows the RISC engine enough free time to issue commands and check the status of the other units in addition to possibly doing some calculations itself to update the data structure or determine the segment to process next.
- the RISC processor issues instructions to the functional units via memory mapped writes.
- the address of the write specifies the functional unit along with the opcode of the instruction for that functional unit.
- the 32-bits of data provided by the write are parameters of the opcode.
- Each opcode specifies how it interprets the data bus.
- Each functional unit can only handle one instruction at a time. If a unit is busy processing an instruction previously issued, any new instructions will be pended at the PCU.
- the RISC processor is also free to read registers and status in each of the functional units to determine the current state of the machine.
- an instruction queue is provided in the PCU to queue up a number of instructions.
- the instruction at the head of the queue is dispatched to the appropriate functional unit as soon as that unit is free to take a new instruction. If the unit is busy, then the queue will not issue the instruction and will stall all subsequent instructions.
- This queuing mechanism allows the dependencies between the function units to be encoded in the ordering of the queued instructions. At any time, the processor is free to by-pass the instruction queue and issue instructions directly to and read information from the function units.
- viscosity of the simulation could be modified using the over-relaxation technique described by Traub et al. in "Viscosity Reduction in Physical Process Simulation", U.S. application Ser. No. 08/255,409, filed Jun. 8, 1994.
- Viscosity is a measure of a fluid's resistance to a shear force (i.e., a force which acts parallel to the direction of fluid flow).
- a shear force i.e., a force which acts parallel to the direction of fluid flow.
- viscosity results from interactions between neighboring particles in the fluid that cause the velocities of the particles to gravitate toward an average value.
- viscosity results from interactions between particles positioned in specific voxels that cause the net velocity of the particles positioned in a voxel to gravitate toward the net velocity of the particles positioned in neighboring voxels.
- each voxel in a lattice system represents a region of simulated space that is substantially larger than the physical space that would be occupied by an actual particle
- the viscosity resulting from interactions between voxels is substantially greater than that resulting from molecular particle interactions in real fluids (i.e., the "averaging" resulting from each voxel interaction affects a substantially larger region of space than that resulting from each molecular particle interaction).
- Viscosity in a lattice system can be reduced by increasing the density of the lattice (i.e., by decreasing the quantity of simulated space that is represented by each voxel), and can also be reduced through use of over-relaxation.
- Viscosity, ⁇ can be expressed in terms of ⁇ , the relaxation parameter: ##EQU27## where T is the temperature of the fluid.
- T the temperature of the fluid.
- over-relaxation effectively increases the density of the lattice.
- Use of over-relaxation therefore has a dramatic effect on the processing necessary to simulate a physical system with a particular resolution (or the resolution with which a particular processor can simulate a physical system). For example, a tenfold increase in the effective density of a three dimensional lattice reduces the processing required to simulate a physical system with the lattice to a particular level of resolution by a factor of almost ten thousand (i.e., ten cubed less the additional processing required to implement over-relaxation and multiplied by a tenfold decrease in the time required to simulate a fluid of a given velocity).
- the computer system performs viscosity modification operations on the state vectors. These operations are typically performed after the microdynamics operations and apply a set of rules that are similar to, or the same as, the rules applied during the microdynamics operations. Where the same rules are applied, the rules modify the state vectors by a first amount during the microdynamics operations and a second amount during the viscosity modification operations, where the first amount is related to the second amount by a relaxation parameter. Because the rules used in the interaction operations conserve mass, momentum and energy, this approach ensures that these properties will be conserved during the viscosity modification operations.
- the result of the multiplication is truncated to ensure that the second amount is an integer value.
- a random value between zero and one is added to the result of the multiplication prior to truncation.
- the viscosity of the lattice system is reduced by using a relaxation parameter having a value greater than one and less than two.
- the relaxation parameter approaches two, the viscosity of the simulated system approaches zero and the system becomes unstable.
- Viscosity which is essentially a form of friction, tends to damp out fluctuations in the system.
- instability occurs when there is no viscosity because these fluctuations are allowed to spread unchecked through the system. It has been found that instability can generally be avoided by using a relaxation parameter that is less than or equal to 1.9.
Abstract
Description
Re=uL/v.
c.sub.i =(c.sub.x, c.sub.y, c.sub.z, c.sub.w)
S={F.sub.α }
V.sub.iα =|c.sub.i ·n.sub.α |A.sub.α.
Γ.sub.iα (x)=N.sub.i (x)V.sub.iα (x).
Γ.sub.iα (x)=N.sub.i (x)V.sub.iα (x).
Γ.sub.iα (β)=Γ.sub.i (β)V.sub.iα (β)/V.sub.iα,
N.sub.i (α)=Γ.sub.iIN (α)/V.sub.iα,
N.sub.i (α)=Γ.sub.iOTHER (α)/V.sub.iα,
Γ.sub.iOTHER (α, t)=Γ.sub.iOUT (α, t-1).
P.sub.n (α)=n.sub.α ·P(α).
ΔΓ.sub.i (α)=Γ.sub.iIN (α)-N.sub.n-Bi (α)V.sub.iα.
Γ.sub.iOUT (α)=N.sub.n-Bi (α)V.sub.iα -ΔΓ.sub.i* (α),
P.sub.n (α)=n.sub.α ·P(α).
u.sub.t (α)=(P(α)-P.sub.n (α)n.sub.α)/ρ,
<N>.sub.j.sup.(eq) (j=0,1,2)
ΔΓ.sub.i (α)=Γ.sub.iIN (α)-N.sub.n-Bi (α)V.sub.iα.
Γ.sub.iOUT (α)=N.sub.n-Bi (α)V.sub.iα -ΔΓ.sub.i* (α)+K.sub.ƒ (n.sub.α ·c.sub.i) N.sub.n-Bi* (α)-N.sub.n-Bi (α)!V.sub.iα,
Γ.sub.αjiOUT.sbsb.f =Γ.sub.αjiOUT +δΓ.sub.αji
N.sub.i (x, t.sub.s +)=floor (scale(N.sub.i (x, t.sub.s -)+N.sub.iF→V (x))+rand)/|scale!
c.sub.i +c.sub.j =c.sub.k +c.sub.1
c.sub.i.sup.2 +c.sub.j.sup.2 =c.sub.k.sup.2 +c.sub.1.sup.2
N.sub.scatt =A.sub.j R.sub.f *N.sub.i *N.sub.j -R.sub.b *N.sub.k *N.sub.l!
<N>.sub.j.sup.(eq) (j=0,1,2)
N.sub.scatt =floor A.sub.js *(N.sub.i *N.sub.j -N.sub.k *N.sub.l)+rand)/scale!
N.sub.scatt =floor (R.sub.ƒ '*N.sub.i *N.sub.j -R.sub.b '*N.sub.k *N.sub.l)+rand)/scale!
R.sub.f *N.sub.i *N.sub.j =R.sub.b *N.sub.k *N.sub.l
ln r.sub.i.sup.-1 N.sub.i +ln r.sub.j.sup.-1 N.sub.j =ln r.sub.k.sup.-1 N.sub.k +ln r.sub.1.sup.-1 N.sub.1
r.sub.i.sup.-1 r.sub.j.sup.-1 =R.sub.f
r.sub.k.sup.-1 r.sub.l.sup.-1 =R.sub.b.
V.sub.iα =|c.sub.i ·n.sub.α |A.sub.α
Claims (35)
Priority Applications (4)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
US08/566,032 US5848260A (en) | 1993-12-10 | 1995-12-01 | Computer system for simulating physical processes |
PCT/US1996/018452 WO1997021195A1 (en) | 1995-12-01 | 1996-11-15 | Computer system for simulating physical processes |
EP96940800A EP0870283A1 (en) | 1995-12-01 | 1996-11-15 | Computer system for simulating physical processes |
JP52128997A JP3932142B2 (en) | 1995-12-01 | 1996-11-15 | Computer system for simulating physical processes |
Applications Claiming Priority (2)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
US08/165,293 US5594671A (en) | 1993-03-12 | 1993-12-10 | Computer system for simulating physical processes using multiple-integer state vectors |
US08/566,032 US5848260A (en) | 1993-12-10 | 1995-12-01 | Computer system for simulating physical processes |
Related Parent Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
US08/165,293 Continuation-In-Part US5594671A (en) | 1993-03-12 | 1993-12-10 | Computer system for simulating physical processes using multiple-integer state vectors |
Publications (1)
Publication Number | Publication Date |
---|---|
US5848260A true US5848260A (en) | 1998-12-08 |
Family
ID=24261178
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
US08/566,032 Expired - Lifetime US5848260A (en) | 1993-12-10 | 1995-12-01 | Computer system for simulating physical processes |
Country Status (4)
Country | Link |
---|---|
US (1) | US5848260A (en) |
EP (1) | EP0870283A1 (en) |
JP (1) | JP3932142B2 (en) |
WO (1) | WO1997021195A1 (en) |
Cited By (19)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2001013256A1 (en) * | 1999-08-12 | 2001-02-22 | Margolus Norman H | A mechanism for efficient data access and communication in parallel computations on an emulated spatial lattice |
US20030039159A1 (en) * | 2001-02-28 | 2003-02-27 | Ware Frederick A. | Asynchronous, high-bandwidth memory component using calibrated timing elements |
US6556961B1 (en) | 1999-12-30 | 2003-04-29 | Innovative Computing Group, Inc. | Method and apparatus for process modeling using a dynamical system |
US20060235562A1 (en) * | 2001-09-26 | 2006-10-19 | Nederlandse Organisatie Voor Toegepast-Natuurwetenschappelijk Onderzoek Tno | Method and system for manufacturing a three-dimensional object |
US20080004838A1 (en) * | 2005-12-08 | 2008-01-03 | The Johns Hopkins University | Air-Speed Wind Tunnel Data Analysis Suite |
US20080109105A1 (en) * | 2005-05-27 | 2008-05-08 | Mueller Weingarten Ag | Method for optimizing the transport displacement of workpieces in transfer presses |
US20090070079A1 (en) * | 2007-09-11 | 2009-03-12 | Prometech Software, Inc. | Sliced data structure for particle-based simulation, and method for loading particle-based simulation using sliced data structure into GPU |
US20100088081A1 (en) * | 2006-08-10 | 2010-04-08 | Exa Corporation | Computer Simulation of Physical Processes |
US7948485B1 (en) * | 2005-12-12 | 2011-05-24 | Sony Computer Entertainment Inc. | Real-time computer simulation of water surfaces |
US20130027407A1 (en) * | 2011-07-27 | 2013-01-31 | Dreamworks Animation Llc | Fluid dynamics framework for animated special effects |
WO2015017648A3 (en) * | 2013-07-31 | 2015-11-05 | Exa Corporation | Temperature coupling algorithm for hybrid thermal lattice boltzmann method |
US20150356217A1 (en) * | 2013-07-24 | 2015-12-10 | Exa Corportion | Lattice Boltzmann Collision Operators Enforcing Isotropy and Galilean Invariance |
CN111428423A (en) * | 2019-01-10 | 2020-07-17 | 达索系统西姆利亚公司 | Lattice boltzmann solver for realizing total energy conservation |
EP3933654A1 (en) | 2020-06-29 | 2022-01-05 | Dassault Systemes Simulia Corp. | Computer system for simulating physical processes using surface algorithm |
US11461512B2 (en) | 2017-01-26 | 2022-10-04 | Dassault Systemes Simulia Corp. | Multi-phase flow visualizations based on fluid occupation time |
US11530598B2 (en) | 2018-08-21 | 2022-12-20 | Dassault Systemes Simulia Corp. | Determination of oil removed by gas via miscible displacement in reservoir rock |
US11613984B2 (en) | 2019-09-04 | 2023-03-28 | Dassault Systemes Simulia Corp. | Determination of hydrocarbon mobilization potential for enhanced oil recovery |
US11714040B2 (en) | 2018-01-10 | 2023-08-01 | Dassault Systemes Simulia Corp. | Determining fluid flow characteristics of porous mediums |
US11907625B2 (en) | 2020-12-29 | 2024-02-20 | Dassault Systemes Americas Corp. | Computer simulation of multi-phase and multi-component fluid flows including physics of under-resolved porous structures |
Families Citing this family (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US5910902A (en) * | 1997-03-28 | 1999-06-08 | Exa Corporation | Computer simulation of physical processes |
US5953239A (en) * | 1997-12-29 | 1999-09-14 | Exa Corporation | Computer simulation of physical processes |
US6089744A (en) * | 1997-12-29 | 2000-07-18 | Exa Corporation | Computer simulation of physical processes |
JP4615543B2 (en) * | 2000-12-12 | 2011-01-19 | 富士通株式会社 | Coupling analysis method and program thereof |
EP1260920A1 (en) * | 2001-05-23 | 2002-11-27 | Fraunhofer-Gesellschaft zur Förderung der angewandten Forschung e.V. | Method and apparatus for computing an interface of a fluid in a space |
US11379636B2 (en) | 2019-01-10 | 2022-07-05 | Dassault Systemes Simulia Corp. | Lattice Boltzmann solver enforcing total energy conservation |
Citations (21)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US4498134A (en) * | 1982-01-26 | 1985-02-05 | Hughes Aircraft Company | Segregator functional plane for use in a modular array processor |
US4507726A (en) * | 1982-01-26 | 1985-03-26 | Hughes Aircraft Company | Array processor architecture utilizing modular elemental processors |
EP0228915A2 (en) * | 1985-12-27 | 1987-07-15 | Thinking Machines Corporation | Method and apparatus for simulating systems described by partial differential equations |
US4729098A (en) * | 1985-06-05 | 1988-03-01 | General Electric Company | System and method employing nonlinear interpolation for the display of surface structures contained within the interior region of a solid body |
US4831519A (en) * | 1985-12-12 | 1989-05-16 | Itt Corporation | Cellular array processor with variable nesting depth vector control by selective enabling of left and right neighboring processor cells |
US4860245A (en) * | 1986-10-08 | 1989-08-22 | Nec Corporation | Vector data processing system comprising an I/O control for each vector data processor and another I/O control for at least one other vector data processor |
US4933895A (en) * | 1987-07-10 | 1990-06-12 | Hughes Aircraft Company | Cellular array having data dependent processing capabilities |
US4969116A (en) * | 1987-09-29 | 1990-11-06 | Kabushiki Kaisha Toshiba | Method of generating discretization grid for finite-difference simulation |
US4989166A (en) * | 1988-07-12 | 1991-01-29 | Hitachi, Ltd. | Method for synthesizing analysis model and flow analysis system |
US5038302A (en) * | 1988-07-26 | 1991-08-06 | The Research Foundation Of State University Of New York | Method of converting continuous three-dimensional geometrical representations into discrete three-dimensional voxel-based representations within a three-dimensional voxel-based system |
WO1992001993A1 (en) * | 1990-07-12 | 1992-02-06 | Massachusetts Institute Of Technology | A fluid dynamics parallel processor |
US5159690A (en) * | 1988-09-30 | 1992-10-27 | Massachusetts Institute Of Technology | Multidimensional cellular data array processing system which separately permutes stored data elements and applies transformation rules to permuted elements |
US5255212A (en) * | 1990-03-27 | 1993-10-19 | Yoshiomi Kondoh | Method of predicting a physical quantity of a fluid or a magnetofluid |
US5361385A (en) * | 1992-08-26 | 1994-11-01 | Reuven Bakalash | Parallel computing system for volumetric modeling, data processing and visualization |
US5377129A (en) * | 1990-07-12 | 1994-12-27 | Massachusetts Institute Of Technology | Particle interaction processing system |
US5408638A (en) * | 1990-12-21 | 1995-04-18 | Hitachi, Ltd. | Method of generating partial differential equations for simulation, simulation method, and method of generating simulation programs |
US5416729A (en) * | 1992-06-24 | 1995-05-16 | Nippon Telegraph And Telephone Corporation | Generalized solids modeling for three-dimensional topography simulation |
US5424963A (en) * | 1992-11-25 | 1995-06-13 | Photon Research Associates, Inc. | Molecular dynamics simulation method and apparatus |
US5442733A (en) * | 1992-03-20 | 1995-08-15 | The Research Foundation Of State University Of New York | Method and apparatus for generating realistic images using a discrete representation |
US5548694A (en) * | 1995-01-31 | 1996-08-20 | Mitsubishi Electric Information Technology Center America, Inc. | Collision avoidance system for voxel-based object representation |
US5606517A (en) * | 1994-06-08 | 1997-02-25 | Exa Corporation | Viscosity reduction in physical process simulation |
-
1995
- 1995-12-01 US US08/566,032 patent/US5848260A/en not_active Expired - Lifetime
-
1996
- 1996-11-15 JP JP52128997A patent/JP3932142B2/en not_active Expired - Lifetime
- 1996-11-15 EP EP96940800A patent/EP0870283A1/en not_active Withdrawn
- 1996-11-15 WO PCT/US1996/018452 patent/WO1997021195A1/en not_active Application Discontinuation
Patent Citations (24)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US4498134A (en) * | 1982-01-26 | 1985-02-05 | Hughes Aircraft Company | Segregator functional plane for use in a modular array processor |
US4507726A (en) * | 1982-01-26 | 1985-03-26 | Hughes Aircraft Company | Array processor architecture utilizing modular elemental processors |
US4729098A (en) * | 1985-06-05 | 1988-03-01 | General Electric Company | System and method employing nonlinear interpolation for the display of surface structures contained within the interior region of a solid body |
US4831519A (en) * | 1985-12-12 | 1989-05-16 | Itt Corporation | Cellular array processor with variable nesting depth vector control by selective enabling of left and right neighboring processor cells |
EP0228915A2 (en) * | 1985-12-27 | 1987-07-15 | Thinking Machines Corporation | Method and apparatus for simulating systems described by partial differential equations |
US4809202A (en) * | 1985-12-27 | 1989-02-28 | Thinking Machines Corporation | Method and apparatus for simulating systems described by partial differential equations |
US4860245A (en) * | 1986-10-08 | 1989-08-22 | Nec Corporation | Vector data processing system comprising an I/O control for each vector data processor and another I/O control for at least one other vector data processor |
US4933895A (en) * | 1987-07-10 | 1990-06-12 | Hughes Aircraft Company | Cellular array having data dependent processing capabilities |
US4969116A (en) * | 1987-09-29 | 1990-11-06 | Kabushiki Kaisha Toshiba | Method of generating discretization grid for finite-difference simulation |
US4989166A (en) * | 1988-07-12 | 1991-01-29 | Hitachi, Ltd. | Method for synthesizing analysis model and flow analysis system |
US5038302A (en) * | 1988-07-26 | 1991-08-06 | The Research Foundation Of State University Of New York | Method of converting continuous three-dimensional geometrical representations into discrete three-dimensional voxel-based representations within a three-dimensional voxel-based system |
US5159690A (en) * | 1988-09-30 | 1992-10-27 | Massachusetts Institute Of Technology | Multidimensional cellular data array processing system which separately permutes stored data elements and applies transformation rules to permuted elements |
US5255212A (en) * | 1990-03-27 | 1993-10-19 | Yoshiomi Kondoh | Method of predicting a physical quantity of a fluid or a magnetofluid |
WO1992001993A1 (en) * | 1990-07-12 | 1992-02-06 | Massachusetts Institute Of Technology | A fluid dynamics parallel processor |
US5377129A (en) * | 1990-07-12 | 1994-12-27 | Massachusetts Institute Of Technology | Particle interaction processing system |
US5432718A (en) * | 1990-07-12 | 1995-07-11 | Massachusetts Institute Of Technology | Particle interaction processing system |
US5408638A (en) * | 1990-12-21 | 1995-04-18 | Hitachi, Ltd. | Method of generating partial differential equations for simulation, simulation method, and method of generating simulation programs |
US5442733A (en) * | 1992-03-20 | 1995-08-15 | The Research Foundation Of State University Of New York | Method and apparatus for generating realistic images using a discrete representation |
US5416729A (en) * | 1992-06-24 | 1995-05-16 | Nippon Telegraph And Telephone Corporation | Generalized solids modeling for three-dimensional topography simulation |
US5361385A (en) * | 1992-08-26 | 1994-11-01 | Reuven Bakalash | Parallel computing system for volumetric modeling, data processing and visualization |
US5424963A (en) * | 1992-11-25 | 1995-06-13 | Photon Research Associates, Inc. | Molecular dynamics simulation method and apparatus |
US5594671A (en) * | 1993-03-12 | 1997-01-14 | Exa Corporation | Computer system for simulating physical processes using multiple-integer state vectors |
US5606517A (en) * | 1994-06-08 | 1997-02-25 | Exa Corporation | Viscosity reduction in physical process simulation |
US5548694A (en) * | 1995-01-31 | 1996-08-20 | Mitsubishi Electric Information Technology Center America, Inc. | Collision avoidance system for voxel-based object representation |
Non-Patent Citations (68)
Title |
---|
Burges et al., "Buoyant Mixtures of Cellular Automation Gases," Complex Systems, 1:31-50 (1987). |
Burges et al., Buoyant Mixtures of Cellular Automation Gases, Complex Systems, 1:31 50 (1987). * |
Chen et al., "Lattice Boltzmann Model for Stimulation of Magnetohydordynamcis," Physical Review Letters, vol. 67, No. 27 (Dec. 1991), pp. 3776-3779. |
Chen et al., Lattice Boltzmann Model for Stimulation of Magnetohydordynamcis, Physical Review Letters, vol. 67, No. 27 (Dec. 1991), pp. 3776 3779. * |
Chopard et al., "Cellular Automata Model for Heat Conduction in a Fluid," Pysics Letters A, vol. 126, Nos. 8, 9, pp. 476-480 (Jan. 1988). |
Chopard et al., Cellular Automata Model for Heat Conduction in a Fluid, Pysics Letters A, vol. 126, Nos. 8, 9, pp. 476 480 (Jan. 1988). * |
Clouqueur et al., "RAPI, A Cellular Automation Machine for Fluid Dynamics," Complex Systems, pp. 585-597 (1987). |
Clouqueur et al., RAPI, A Cellular Automation Machine for Fluid Dynamics, Complex Systems, pp. 585 597 (1987). * |
d Humieres et al., 2 D and 3 D Hydrodynamics on Lattice Gases, Fruhjahrstagung der Schweiz. Physikalischen Gesellschaft, vol. 59 (1986). * |
d Humieres et al., Lattice Gas Models for 3D Hydrodynamics, 2 (4), pp. 291 297 (1986). * |
d Humieres et al., Numerial Simulations of Hydrodynamics with Lattice Gas Automata in Two Dimensions, Complex Systems 1:599 632 (1987). * |
d'Humieres et al., "2-D and 3-D Hydrodynamics on Lattice Gases," Fruhjahrstagung der Schweiz. Physikalischen Gesellschaft, vol. 59 (1986). |
d'Humieres et al., "Lattice Gas Models for 3D Hydrodynamics," 2 (4), pp. 291-297 (1986). |
d'Humieres et al., "Numerial Simulations of Hydrodynamics with Lattice Gas Automata in Two Dimensions," Complex Systems 1:599-632 (1987). |
Fredkin et al., "Conservative Logic," International Journal of Theoretical Pysics, vol. 21, Nos. 3/4, pp. 219-253 (1982). |
Fredkin et al., Conservative Logic, International Journal of Theoretical Pysics, vol. 21, Nos. 3/4, pp. 219 253 (1982). * |
Frisch et al., "Lattice Gas Hydrodynamics in Two and Three Dimensions," Complex Systems 1:649-707 (1987). |
Frisch et al., "Lattice-Gas Automata for the Navier-Stokes Equation," Physical Review Letters, vol. 56, No. 14, pp. 1505-1508 (Apr. 1986). |
Frisch et al., Lattice Gas Automata for the Navier Stokes Equation, Physical Review Letters, vol. 56, No. 14, pp. 1505 1508 (Apr. 1986). * |
Frisch et al., Lattice Gas Hydrodynamics in Two and Three Dimensions, Complex Systems 1:649 707 (1987). * |
Hardy et al., "Molecular Dynamics of a Classical Lattice Gas: Transport Properties and Time Correlation Functions," Physical Review A, vol.13, No. 5, pp. 1949-1961 (May 1976). |
Hardy et al., Molecular Dynamics of a Classical Lattice Gas: Transport Properties and Time Correlation Functions, Physical Review A, vol.13, No. 5, pp. 1949 1961 (May 1976). * |
Hasslacher, "Discrete Fluids: Part I: Background for Lattice Gas Automata," Los Alamos Science Special Issue, pp. 175-217 (1987). |
Hasslacher, Discrete Fluids: Part I: Background for Lattice Gas Automata, Los Alamos Science Special Issue, pp. 175 217 (1987). * |
Henon, "Isometric Collision Rules for the Four-Dimensional FCHC Lattice Gas," Complex Systems 1:475-494 (1987). |
Henon, Isometric Collision Rules for the Four Dimensional FCHC Lattice Gas, Complex Systems 1:475 494 (1987). * |
Kadanoff et al., "Transport Coefficients Near the Critical Point: A Master Equation Approach," 165 Phys. Review 310 (1968). |
Kadanoff et al., Transport Coefficients Near the Critical Point: A Master Equation Approach, 165 Phys. Review 310 (1968). * |
Kadanoff, "On Two Levels" Physics Today (Sep. 1986). |
Kadanoff, On Two Levels Physics Today (Sep. 1986). * |
Kaufman et al "Volume Graphics," IEEE, Jul. 1993, pp. 51-64. |
Kaufman et al Volume Graphics, IEEE, Jul. 1993, pp. 51 64. * |
Klein, "Negative Absolute Temperatures," Physical Review, vol. 104, No. 3, p. 589 (Nov. 1956). |
Klein, Negative Absolute Temperatures, Physical Review, vol. 104, No. 3, p. 589 (Nov. 1956). * |
Lawerence J. Rosenblum, "Photorealistic Terrain Imaging and Flight Simulation," IEEE, 1994 pp. 10-12. |
Lawerence J. Rosenblum, Photorealistic Terrain Imaging and Flight Simulation, IEEE, 1994 pp. 10 12. * |
Lin et al., "Plasma Simulation Using the Massively Parallel Processor," NASA CP-2478, pp. 185-191 (Jul. 1987). |
Lin et al., Plasma Simulation Using the Massively Parallel Processor, NASA CP 2478, pp. 185 191 (Jul. 1987). * |
Margolus et al., "Cellular-Automata Supercomputers for Fluid-Dynamics Modeling," Physical Review Letters, vol. 56, No. 16, pp. 1694-1696 (Apr. 1986). |
Margolus et al., Cellular Automata Supercomputers for Fluid Dynamics Modeling, Physical Review Letters, vol. 56, No. 16, pp. 1694 1696 (Apr. 1986). * |
Molvig et al., "Multi-species Lattice-Gas Automata for Realistic Fluid Dynamics," Springer Proceedings in Pysics, vol. 46, pp. 206-231 (1990). |
Molvig et al., "Removing the Discreteness Artifacts in 3D Lattice-Gas Fluids," Proceedings of the Workshop on Discrete Kinetic Theory. Lattice Gas Dynamics and Foundations of Hydrodynamics (1989). |
Molvig et al., Multi species Lattice Gas Automata for Realistic Fluid Dynamics, Springer Proceedings in Pysics, vol. 46, pp. 206 231 (1990). * |
Molvig et al., Removing the Discreteness Artifacts in 3D Lattice Gas Fluids, Proceedings of the Workshop on Discrete Kinetic Theory. Lattice Gas Dynamics and Foundations of Hydrodynamics (1989). * |
Perez et al, "OUPPI-1, A SIMD Computer Using Integrated Parallel Processors," Proc. of Conpar Conf. Papers Plenary Sess. and Stream A Manchester, Sep. 12-16, 1988, British Comp. Soc. Parallel Proc. Spec. Group GB. |
Perez et al, OUPPI 1, A SIMD Computer Using Integrated Parallel Processors, Proc. of Conpar Conf. Papers Plenary Sess. and Stream A Manchester, Sep. 12 16, 1988, British Comp. Soc. Parallel Proc. Spec. Group GB. * |
Product Literature, "Digital Physics Technology Principles & Applications," Exa Corporation (Aug. 1994). |
Product Literature, "Understanding Fluid CAD, ExaResolute Product Technical Description," Exa Corporation (1994). |
Product Literature, Digital Physics Technology Principles & Applications, Exa Corporation (Aug. 1994). * |
Product Literature, Understanding Fluid CAD, ExaResolute Product Technical Description, Exa Corporation (1994). * |
Qian et al., "Lattice BGK Models for Navier-Stokes Equation," Europhysics Letters, 17(6):479-484 (1992). |
Qian et al., Lattice BGK Models for Navier Stokes Equation, Europhysics Letters, 17(6):479 484 (1992). * |
Ramsey, "Thermodynamics and Statistical Mechanics at Negative Absolute Temperatures," Physical Review, vol. 103, No. 1, pp. 20-28 (Jul. 1956). |
Ramsey, Thermodynamics and Statistical Mechanics at Negative Absolute Temperatures, Physical Review, vol. 103, No. 1, pp. 20 28 (Jul. 1956). * |
Rivet, "Three-Dimensional Lattice Gas Hydrodynamical Simulations: First Results," C.R. Acad. Sci. Paris, t. 305, Series II, pp. 751-756 (1987). |
Rivet, Three Dimensional Lattice Gas Hydrodynamical Simulations: First Results, C.R. Acad. Sci. Paris, t. 305, Series II, pp. 751 756 (1987). * |
Teixeira, "Continuum Limit of Lattice Gas Fluid Dynamics," B.A. Sc. Engineering Science University of Toronto (1988). |
Teixeira, Continuum Limit of Lattice Gas Fluid Dynamics, B.A. Sc. Engineering Science University of Toronto (1988). * |
Toffoli, "Cellular Automata as an Alternative to (Rather Than an Approximation of) Differentiation Equation in Modeling Physics," Physica 10D:117-127 (1984). |
Toffoli, Cellular Automata as an Alternative to (Rather Than an Approximation of) Differentiation Equation in Modeling Physics, Physica 10D:117 127 (1984). * |
Tucker, "Cellular Automata Machine: The Ultimate Parallel Computer," High Technology, pp. 85-87 (Jun. 1984). |
Tucker, Cellular Automata Machine: The Ultimate Parallel Computer, High Technology, pp. 85 87 (Jun. 1984). * |
Vichniac, "Cellular-Automata Fluids," Instabilities and Nonequilibrium Structures II, pp. 97-116 (1989). |
Vichniac, "Stimulating Physics with Cellular Automata," Physica, 10D:96-110 (1984). |
Vichniac, Cellular Automata Fluids, Instabilities and Nonequilibrium Structures II, pp. 97 116 (1989). * |
Vichniac, Stimulating Physics with Cellular Automata, Physica, 10D:96 110 (1984). * |
Wolfram, "Cellular Automation Fluids 1: Basic Theory," J. Stat. Phys., vol. 45, Nos. 3/4, pp. 471-526 (1986). |
Wolfram, Cellular Automation Fluids 1: Basic Theory, J. Stat. Phys., vol. 45, Nos. 3/4, pp. 471 526 (1986). * |
Cited By (45)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US6205533B1 (en) * | 1999-08-12 | 2001-03-20 | Norman H. Margolus | Mechanism for efficient data access and communication in parallel computations on an emulated spatial lattice |
WO2001013256A1 (en) * | 1999-08-12 | 2001-02-22 | Margolus Norman H | A mechanism for efficient data access and communication in parallel computations on an emulated spatial lattice |
US6556961B1 (en) | 1999-12-30 | 2003-04-29 | Innovative Computing Group, Inc. | Method and apparatus for process modeling using a dynamical system |
US8295107B2 (en) | 2001-02-28 | 2012-10-23 | Rambus Inc. | Asynchronous pipelined memory access |
US7362626B2 (en) | 2001-02-28 | 2008-04-22 | Rambus Inc. | Asynchronous, high-bandwidth memory component using calibrated timing elements |
US6788593B2 (en) * | 2001-02-28 | 2004-09-07 | Rambus, Inc. | Asynchronous, high-bandwidth memory component using calibrated timing elements |
US6788594B2 (en) | 2001-02-28 | 2004-09-07 | Rambus Inc. | Asynchronous, high-bandwidth memory component using calibrated timing elements |
US20040213052A1 (en) * | 2001-02-28 | 2004-10-28 | Ware Frederick A | Asynchronous, high-bandwidth memory component using calibrated timing elements |
US6934201B2 (en) | 2001-02-28 | 2005-08-23 | Rambus Inc. | Asynchronous, high-bandwidth memory component using calibrated timing elements |
US20050237851A1 (en) * | 2001-02-28 | 2005-10-27 | Ware Frederick A | Asynchronous, high-bandwidth memory component using calibrated timing elements |
US7830735B2 (en) | 2001-02-28 | 2010-11-09 | Rambus Inc. | Asynchronous, high-bandwidth memory component using calibrated timing elements |
US8059476B2 (en) | 2001-02-28 | 2011-11-15 | Rambus Inc. | Control component for controlling a delay interval within a memory component |
US6574153B2 (en) | 2001-02-28 | 2003-06-03 | Rambus Inc. | Asynchronous, high-bandwidth memory component using calibrated timing elements |
US20090213670A1 (en) * | 2001-02-28 | 2009-08-27 | Ware Frederick A | Asynchronous, high-bandwidth memory component using calibrated timing elements |
US20080144408A1 (en) * | 2001-02-28 | 2008-06-19 | Ware Frederick A | Asynchronous, high-bandwidth memory component using calibrated timing elements |
US20030039159A1 (en) * | 2001-02-28 | 2003-02-27 | Ware Frederick A. | Asynchronous, high-bandwidth memory component using calibrated timing elements |
US7529141B2 (en) | 2001-02-28 | 2009-05-05 | Rambus Inc. | Asynchronous, high-bandwidth memory component using calibrated timing elements |
US20060235562A1 (en) * | 2001-09-26 | 2006-10-19 | Nederlandse Organisatie Voor Toegepast-Natuurwetenschappelijk Onderzoek Tno | Method and system for manufacturing a three-dimensional object |
US20080109105A1 (en) * | 2005-05-27 | 2008-05-08 | Mueller Weingarten Ag | Method for optimizing the transport displacement of workpieces in transfer presses |
US7636612B2 (en) * | 2005-05-27 | 2009-12-22 | Mueller Weingarten Ag | Method for optimizing the transport displacement of workpieces in transfer presses |
US20080004838A1 (en) * | 2005-12-08 | 2008-01-03 | The Johns Hopkins University | Air-Speed Wind Tunnel Data Analysis Suite |
US7948485B1 (en) * | 2005-12-12 | 2011-05-24 | Sony Computer Entertainment Inc. | Real-time computer simulation of water surfaces |
US20100088081A1 (en) * | 2006-08-10 | 2010-04-08 | Exa Corporation | Computer Simulation of Physical Processes |
US8224633B2 (en) * | 2006-08-10 | 2012-07-17 | Exa Corporation | Computer simulation of physical processes |
US8346522B2 (en) | 2006-08-10 | 2013-01-01 | Exa Corporation | Computer simulation of physical processes |
US7920996B2 (en) * | 2007-09-11 | 2011-04-05 | Prometech Software, Inc. | Sliced data structure for particle-based simulation, and method for loading particle-based simulation using sliced data structure into GPU |
US20090070079A1 (en) * | 2007-09-11 | 2009-03-12 | Prometech Software, Inc. | Sliced data structure for particle-based simulation, and method for loading particle-based simulation using sliced data structure into GPU |
US9984489B2 (en) * | 2011-07-27 | 2018-05-29 | Dreamworks Animation L.L.C. | Fluid dynamics framework for animated special effects |
US20130027407A1 (en) * | 2011-07-27 | 2013-01-31 | Dreamworks Animation Llc | Fluid dynamics framework for animated special effects |
US11194941B2 (en) * | 2013-07-24 | 2021-12-07 | Dassault Systemes Simulia Corp. | Lattice Boltzmann collision operators enforcing isotropy and Galilean invariance |
US9576087B2 (en) * | 2013-07-24 | 2017-02-21 | Exa Corporation | Lattice Boltzmann collision operators enforcing isotropy and Galilean invariance |
US20150356217A1 (en) * | 2013-07-24 | 2015-12-10 | Exa Corportion | Lattice Boltzmann Collision Operators Enforcing Isotropy and Galilean Invariance |
US10867088B2 (en) * | 2013-07-24 | 2020-12-15 | Dassault Systemes Simulia Corp. | Lattice boltzmann collision operators enforcing isotropy and galilean invariance |
WO2015017648A3 (en) * | 2013-07-31 | 2015-11-05 | Exa Corporation | Temperature coupling algorithm for hybrid thermal lattice boltzmann method |
US10762252B2 (en) | 2013-07-31 | 2020-09-01 | Dassault Systemes Simulia Corp. | Temperature coupling algorithm for hybrid thermal lattice boltzmann method |
JP2016528628A (en) * | 2013-07-31 | 2016-09-15 | エクサ コーポレイション | Temperature coupling algorithm for hybrid hot-lattice Boltzmann method |
US11941331B2 (en) | 2017-01-26 | 2024-03-26 | Dassault Systemes Americas Corp. | Multi-phase flow visualizations based on fluid occupation time |
US11461512B2 (en) | 2017-01-26 | 2022-10-04 | Dassault Systemes Simulia Corp. | Multi-phase flow visualizations based on fluid occupation time |
US11714040B2 (en) | 2018-01-10 | 2023-08-01 | Dassault Systemes Simulia Corp. | Determining fluid flow characteristics of porous mediums |
US11530598B2 (en) | 2018-08-21 | 2022-12-20 | Dassault Systemes Simulia Corp. | Determination of oil removed by gas via miscible displacement in reservoir rock |
CN111428423A (en) * | 2019-01-10 | 2020-07-17 | 达索系统西姆利亚公司 | Lattice boltzmann solver for realizing total energy conservation |
US11613984B2 (en) | 2019-09-04 | 2023-03-28 | Dassault Systemes Simulia Corp. | Determination of hydrocarbon mobilization potential for enhanced oil recovery |
US11847391B2 (en) | 2020-06-29 | 2023-12-19 | Dassault Systemes Simulia Corp. | Computer system for simulating physical processes using surface algorithm |
EP3933654A1 (en) | 2020-06-29 | 2022-01-05 | Dassault Systemes Simulia Corp. | Computer system for simulating physical processes using surface algorithm |
US11907625B2 (en) | 2020-12-29 | 2024-02-20 | Dassault Systemes Americas Corp. | Computer simulation of multi-phase and multi-component fluid flows including physics of under-resolved porous structures |
Also Published As
Publication number | Publication date |
---|---|
WO1997021195A1 (en) | 1997-06-12 |
EP0870283A4 (en) | 1998-10-14 |
JP3932142B2 (en) | 2007-06-20 |
EP0870283A1 (en) | 1998-10-14 |
JP2000502200A (en) | 2000-02-22 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
US5848260A (en) | Computer system for simulating physical processes | |
US5606517A (en) | Viscosity reduction in physical process simulation | |
US5594671A (en) | Computer system for simulating physical processes using multiple-integer state vectors | |
US5953239A (en) | Computer simulation of physical processes | |
US8346522B2 (en) | Computer simulation of physical processes | |
US10762252B2 (en) | Temperature coupling algorithm for hybrid thermal lattice boltzmann method | |
Mavriplis | Multigrid strategies for viscous flow solvers on anisotropic unstructured meshes | |
US6089744A (en) | Computer simulation of physical processes | |
Mavriplis | Directional agglomeration multigrid techniques for high-Reynolds-number viscous flows | |
US5640335A (en) | Collision operators in physical process simulation | |
US5910902A (en) | Computer simulation of physical processes | |
US10867088B2 (en) | Lattice boltzmann collision operators enforcing isotropy and galilean invariance | |
US11847391B2 (en) | Computer system for simulating physical processes using surface algorithm | |
US11763048B2 (en) | Computer simulation of physical fluids on a mesh in an arbitrary coordinate system | |
Belov | A new implicit multigrid-driven algorithm for unsteady incompressible flow calculations on parallel computers | |
CN112069742A (en) | Stabilizing explicit numerical schemes | |
Ujaldon et al. | The GPU on irregular computing: Performance issues and contributions | |
Flannery | A hybrid fluid simulation on the Graphics Processing Unit (GPU) | |
Strating et al. | Fast computation of stationary inviscid flow around an airfoil | |
Kordulla | On the efficient use of large data bases in the numerical solution of the Navier-Stokes equations on a CRAY computer | |
Madill | Targeted smoke simulation combining control and turbulent flow | |
Bumbacher et al. | A Smoothed Particly Hydrodynamics Approach to the Interactive Swarm Orchestra | |
Umenhoffer et al. | Simulation methods for elastic and fluid materials |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
AS | Assignment |
Owner name: EXA CORPORATION, MASSACHUSETTS Free format text: ASSIGNMENT OF ASSIGNORS INTEREST;ASSIGNORS:CHEN, HUDONG;HILL, DAVID LAWRENCE;HOCH, JAMES;AND OTHERS;REEL/FRAME:007965/0010 Effective date: 19960212 |
|
STCF | Information on status: patent grant |
Free format text: PATENTED CASE |
|
CC | Certificate of correction | ||
FEPP | Fee payment procedure |
Free format text: PAYOR NUMBER ASSIGNED (ORIGINAL EVENT CODE: ASPN); ENTITY STATUS OF PATENT OWNER: SMALL ENTITY |
|
FPAY | Fee payment |
Year of fee payment: 4 |
|
REMI | Maintenance fee reminder mailed | ||
FEPP | Fee payment procedure |
Free format text: PAYOR NUMBER ASSIGNED (ORIGINAL EVENT CODE: ASPN); ENTITY STATUS OF PATENT OWNER: SMALL ENTITY Free format text: PAYER NUMBER DE-ASSIGNED (ORIGINAL EVENT CODE: RMPN); ENTITY STATUS OF PATENT OWNER: SMALL ENTITY |
|
FPAY | Fee payment |
Year of fee payment: 8 |
|
AS | Assignment |
Owner name: SILICON VALLEY BANK,MASSACHUSETTS Free format text: SECURITY AGREEMENT;ASSIGNOR:EXA CORPORATION;REEL/FRAME:024474/0340 Effective date: 20100524 Owner name: SILICON VALLEY BANK, MASSACHUSETTS Free format text: SECURITY AGREEMENT;ASSIGNOR:EXA CORPORATION;REEL/FRAME:024474/0340 Effective date: 20100524 |
|
FPAY | Fee payment |
Year of fee payment: 12 |
|
AS | Assignment |
Owner name: DASSAULT SYSTEMES SIMULIA CORP,, RHODE ISLAND Free format text: MERGER;ASSIGNOR:EXA CORPORATION;REEL/FRAME:048877/0769 Effective date: 20190101 |