HDK
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
GAS_Integrator.h
Go to the documentation of this file.
1 /*
2  * PROPRIETARY INFORMATION. This software is proprietary to
3  * Side Effects Software Inc., and is not to be reproduced,
4  * transmitted, or disclosed in any way without written permission.
5  *
6  * NAME: GAS_Integrator.h ( GAS Library, C++)
7  *
8  * COMMENTS: Used to integrate the forces in a particle fluid
9  * simulation, and to provide several integration
10  * methods.
11  */
12 
13 #ifndef __GAS_Integrator__
14 #define __GAS_Integrator__
15 
16 #include "GAS_API.h"
17 
18 #include "GAS_SubSolver.h"
19 #include "GAS_Utils.h"
20 
21 #include <UT/UT_ValArray.h>
22 #include <UT/UT_Quaternion.h>
23 #include <UT/UT_Vector.h>
24 #include <UT/UT_VectorTypes.h>
26 
27 class GU_Detail;
28 
29 class SIM_GeometryCopy;
30 class SIM_Impacts;
31 class GAS_SPH;
32 
33 // Definitions for number of steps in each integration type
34 #define GAS_RK_45_NUMSTEPS 6
35 #define GAS_RK_4_NUMSTEPS 4
36 #define GAS_RK_21_NUMSTEPS 3
37 #define GAS_RK_32_NUMSTEPS 4
38 #define GAS_EULER_NUMSTEPS 1
39 #define GAS_MIDPOINT_NUMSTEPS 2
40 
42 {
43 public:
45  GET_DATA_FUNC_S(GAS_NAME_GEOMETRY, GeometryName);
46  GET_DATA_FUNC_I("integratetype", IntegrateType);
47  GET_DATA_FUNC_I("advecttype", AdvectType);
48  GET_DATA_FUNC_F("xsphconstant", XsphConstant);
49  GET_DATA_FUNC_B("docollision", DoCollision);
50  GET_DATA_FUNC_B("resolvesdf", ResolveSDF);
51  GET_DATA_FUNC_B("prestashpositions", PreStashPositions);
52  GET_DATA_FUNC_B("collisionuv", CollisionUV);
53  GET_DATA_FUNC_B("docollisionfeedback", DoCollisionFeedback);
54  GET_DATA_FUNC_B("addimpacts", AddImpacts);
55  GET_DATA_FUNC_I("scaletype", ScaleType);
56  GET_DATA_FUNC_B("doincremental", IncrementalForce);
57  GET_DATA_FUNC_B("integrateorientation", IntegrateOrientation);
58  GET_DATA_FUNC_B("integratevelocity", IntegrateVelocity);
59  GET_DATA_FUNC_F("errortolerance", ErrorTolerance);
60  GET_DATA_FUNC_F("repetitiontolerance", RepTolerance);
61  GET_DATA_FUNC_F("minsubstep", MinSubstep);
62  GET_DATA_FUNC_F("maxsubstep", MaxSubstep);
63  GET_DATA_FUNC_B("usepscaleforsdf", UsePScaleForSDF);
64  GET_DATA_FUNC_B("pscaleisradius", PScaleIsRadius);
67 
68  // Available integration types
73 #if 0
74  NEW_RUNGE_KUTTA,
75 #endif
79  LEAPFROG
80  };
81 
82  // Different ways of scaling the forces
83  enum ScaleMethod {
85  DENSITY
86  };
87 
88  // Different advection methods
90  STANDARD = 0,
91  XSPH
92  };
93 
94  enum ObjectState {
95  OBJECTSTATE_INITIAL = 0,
97  OBJECTSTATE_TOTAL
98  };
99 
100  /// Asks all inputs to compute their forces and integrates
101  /// them for the given timestep. We may do this several times
102  /// if a higher-order timestep is used.
103  virtual SIM_Result solveObjectsSubclass(SIM_Engine &engine,
104  SIM_ObjectArray &objects,
105  SIM_ObjectArray &newobjects,
106  SIM_ObjectArray &feedbacktoobjects,
107  const SIM_Time &timestep);
108 
109  bool shouldMultiThread() { return true; }
110 
111 protected:
112  explicit GAS_Integrator(const SIM_DataFactory *factory);
113  virtual ~GAS_Integrator();
114 
115  virtual bool solveGasSubclass(SIM_Engine &engine,
116  SIM_Object *obj,
117  SIM_Time time,
118  SIM_Time timestep)
119  {
120  return true;
121  }
122 
123  virtual void getImpulseMassMatrixSubclass(
124  const SIM_Object &object,
125  const UT_Vector3 &impulseworldpos,
126  UT_DMatrix3 &immatrix) const;
128  const SIM_Object &object,
129  int ptnum, UT_DMatrix3 &immatrix) const;
131  const SIM_Object &object) const;
132 
133 private:
134  /// A struct used to store intermediate information about an object
135  /// during multi-step integration methods.
136  struct ObjectInfo {
137  UT_Vector3Array myInitialPositions;
138  UT_Vector3Array myInitialVelocities;
139  UT_ValArray<UT_Quaternion> myInitialOrientations;
140  UT_Vector3Array myInitialAngularVelocities;
141 
142  UT_Vector3Array myCurrentPositions;
143  UT_Vector3Array myCurrentVelocities;
144  UT_ValArray<UT_Quaternion> myCurrentOrientations;
145  UT_Vector3Array myCurrentAngularVelocities;
146 
147  UT_Vector3Array myTotalPositions;
148  UT_Vector3Array myTotalVelocities;
149  UT_ValArray<UT_Quaternion> myTotalOrientations;
150  UT_Vector3Array myTotalAngularVelocities;
151  };
152 
153  typedef UT_ValArray<UT_Quaternion> gas_QuatArray;
154 
155  /// This data structure is used to store all intermediate information
156  /// used by Runge-Kutta style integrators
157  struct RungeKuttaDataTable {
158  /// We store two arrays of vector arrays - one for positions and
159  /// the other for velocities. Each position/velocity array pair
160  /// represents one step in the particular RK method. Similarly,
161  /// we store orientations/angular velocities as well.
162  UT_Array<UT_Vector3Array> myPositionStepData;
163  UT_Array<UT_Vector3Array> myVelocityStepData;
164  UT_Array<UT_Vector3Array> myOrientationStepData;
165  UT_Array<UT_Vector3Array> myAngVelocityStepData;
166  int myNumSteps;
167 
168  RungeKuttaDataTable() : myNumSteps(0) {}
169 
170  /// Accessor function for each of the step arrays
171  UT_Vector3Array &getPositionStep(int i)
172  {
173  return myPositionStepData(i);
174  }
175 
176  UT_Vector3Array &getVelocityStep(int i)
177  {
178  return myVelocityStepData(i);
179  }
180 
181  UT_Vector3Array &getOrientationStep(int i)
182  {
183  return myOrientationStepData(i);
184  }
185 
186  UT_Vector3Array &getAngVelocityStep(int i)
187  {
188  return myAngVelocityStepData(i);
189  }
190 
191  void resizeIfNeeded(int npts, bool doangular = false)
192  {
193  for (int i = 0; i < myNumSteps; i++)
194  {
195  myPositionStepData(i).entries(npts);
196  myVelocityStepData(i).entries(npts);
197  // Only set these array sizes if doangular is true
198  if (doangular)
199  {
200  myOrientationStepData(i).entries(npts);
201  myAngVelocityStepData(i).entries(npts);
202  }
203  }
204  }
205 
206  void addStepsIfNeeded(int numsteps)
207  {
208  while (myNumSteps < numsteps)
209  {
210  myPositionStepData.append();
211  myVelocityStepData.append();
212  myOrientationStepData.append();
213  myAngVelocityStepData.append();
214 
215  myNumSteps++;
216  }
217  }
218  };
219 
220  /// Solves the intermediate stages for the specified runge-kutta
221  /// integration scheme for each object in 'objects'. The number
222  /// of integration steps is specified in numsteps. 'coefficients'
223  /// is a matrix specifying the runge-kutta coefficients used in
224  /// each integration step. The dimension of coefficients is
225  /// (numsteps-1) x (numsteps-1) since the first step in any RK
226  /// scheme is always the same. 'timecoefficients' are the
227  /// coefficients necessary for incrementing the time at each step.
228  /// This array also has only (numsteps-1) entries, since the first
229  /// of these entries would always be 0 if it had more.
230  void solveRungeKuttaSteps(SIM_Engine &engine,
231  SIM_ObjectArray &objects,
232  SIM_ObjectArray &newobjects,
233  SIM_ObjectArray &feedbacktoobjects,
234  const SIM_Time &timestep,
235  int numsteps,
236  const fpreal64 *coefficients,
237  const fpreal64 *timecoefficients,
238  bool cplUpdate=true);
239 
240  /// This helper function is used to reinitialize each object's
241  /// particle system after a runge-kutta integration step. It
242  /// determines the change from the last step and stores this
243  /// in the appropriate runge-kutta data table. It then sets
244  /// up the system in preparation for the next step.
245  /// 'coefficients' stores the RK coefficients needed to set
246  /// up the state for the next step.
247  void reinitializeSystems(SIM_ObjectArray &objects,
248  int numsteps, int stepnum,
249  const SIM_Time &timestep,
250  const fpreal64 *coefficients,
251  bool cplUpdate=true);
252 
253  class rkInfoParms
254  {
255  public:
256  int numsteps, stepnum;
257  SIM_Time timestep;
258  const fpreal64 *coefficients;
259  bool cplUpdate;
260  };
261  /// Called by the above function to reinitialize a specific
262  /// object's particle system.
263  THREADED_METHOD6(GAS_Integrator, gdp->getNumPoints() > 2048,
264  reinitializeSystem,
265  GU_Detail *, gdp,
266  const GA_PointGroup *, ptgrp,
267  RungeKuttaDataTable &, data,
268  ObjectInfo &, inf,
269  const rkInfoParms &, parms,
270  int, npts)
271  void reinitializeSystemPartial(GU_Detail *gdp,
272  const GA_PointGroup *ptgrp,
273  RungeKuttaDataTable &data,
274  ObjectInfo &inf,
275  const rkInfoParms &parms,
276  int npts,
277  const UT_JobInfo &info);
278 
279  /// This helper function is used to compute the final positions,
280  /// velocities, etc. of all particles using the already-computed
281  /// runge-kutta step data.
282  void computeFinalStates(SIM_ObjectArray &objects,
283  int numsteps, const fpreal64 *coefficients);
284 
285  /// Called by the above function to compute the final state for
286  /// a specific object's particle system.
287  THREADED_METHOD7(GAS_Integrator, gdp->getNumPoints() > 2048,
288  computeFinalState,
289  GU_Detail *, gdp,
290  const GA_PointGroup *, ptgrp,
291  RungeKuttaDataTable &, data,
292  ObjectInfo &, inf,
293  int, numsteps,
294  const fpreal64 *, coefficients,
295  int, npts)
296 
297  void computeFinalStatePartial(GU_Detail *gdp,
298  const GA_PointGroup *ptgrp,
299  RungeKuttaDataTable &data,
300  ObjectInfo &inf, int numsteps,
301  const fpreal64 *coefficients,
302  int npts,
303  const UT_JobInfo &info);
304 
305  /// This is a helper function used to reset everything back to its initial
306  /// state.
307  void resetInitialStates(SIM_ObjectArray &objects);
308 
309  /// Copies the ObjectInfo's current, total, or initial data
310  /// into the geometry data
311  void objectStateToGeo(GU_Detail *gdp,
312  const GA_PointGroup *ptgrp,
313  ObjectInfo &inf,
314  ObjectState state);
315  THREADED_METHOD5(GAS_Integrator, gdp->getNumPoints() > 2048,
316  objectStateToGeo,
317  GU_Detail *, gdp,
318  const GA_PointGroup *, ptgrp,
319  ObjectInfo &, inf,
320  ObjectState, state,
321  GA_Size, npts)
322  void objectStateToGeoPartial(GU_Detail *gdp,
323  const GA_PointGroup *ptgrp,
324  ObjectInfo &inf,
325  ObjectState state,
326  GA_Size npts,
327  const UT_JobInfo &info);
328 
329  /// Returns the maximum number of points to safely iterate over.
330  /// Adds a warning if the point count chaned.
331  GA_Size determinePointCount(const GU_Detail *gdp, ObjectInfo &inf);
332 
333  /// This function just caches the states of all systems in their
334  /// respective object info arrays.
335  void cacheSystemStates(SIM_ObjectArray &objects);
336 
337  /// Copies the gdp's data into the given object state arrays.
338  void geoToObjectState(const GU_Detail *gdp,
339  const GA_PointGroup *ptgrp,
340  ObjectInfo &inf,
341  ObjectState state);
342 
343  THREADED_METHOD5(GAS_Integrator, gdp->getNumPoints() > 2048,
344  geoToObjectState,
345  const GU_Detail *, gdp,
346  const GA_PointGroup *, ptgrp,
347  ObjectInfo &, inf,
348  ObjectState, state,
349  GA_Size, npts)
350  void geoToObjectStatePartial(const GU_Detail *gdp,
351  const GA_PointGroup *ptgrp,
352  ObjectInfo &inf,
353  ObjectState state,
354  GA_Size npts,
355  const UT_JobInfo &info);
356 
357  /// Determines the "distance" between the velociteis currently stored
358  /// in the object geometries and the velocities stored in our object
359  /// info array. We just take the maximum distance for any system.
360  fpreal64 getMaxDistance(SIM_ObjectArray &objects);
361 
362  /// Finds the "distance" described above for the given system.
363  fpreal64 getDistance(const GU_Detail *gdp,
364  const GA_PointGroup *ptgrp,
365  ObjectInfo &inf);
366 
367  /// This function is used for embedded methods to determine whether
368  /// or not an iteration succeeded (ie. if the two
369  /// approximations made by the method agree to a specified accuracy).
370  /// This is done by comparing the myCurrentVelocities array in
371  /// the system's ObjectInfo struct to the current velocity of
372  /// each particle in the system.
373  bool embeddedRKsucceeded(SIM_ObjectArray &objects,
374  const SIM_Time &timestep,
375  fpreal64 toleranceCoefficient);
376 
377  /// Called by the above function on an individual system. This will
378  /// also place the timestep recommended by the embedded method in to
379  /// an attribute on the gdp so that it can be used by GAS_SubStep.
380  bool embeddedRKsucceeded(GU_Detail *gdp,
381  const GA_PointGroup *ptgrp,
382  ObjectInfo &inf,
383  const SIM_Time &timestep, fpreal tolerance,
384  fpreal reptolerance, fpreal minsubstep,
385  fpreal maxsubstep,
386  fpreal64 toleranceCoefficient);
387 
388  /// This function is used for standard integrators to determine
389  /// whether or not an iteration succeeded (ie. if the velocity
390  /// difference (acceleration) does not exceed an error condition)
391  bool integrationSucceeded(SIM_ObjectArray &objects,
392  const SIM_Time &timestep);
393 
394  /// Called by the above function on an individual system. This
395  /// will also place the recommended timestep in to an attribute
396  /// on the gdp so that it can be used by GAS_SubStep
397  bool integrationSucceeded(GU_Detail *gdp,
398  const GA_PointGroup *ptgrp,
399  ObjectInfo &inf,
400  const SIM_Time &timestep, fpreal tolerance,
401  fpreal reptolerance, fpreal minsubstep,
402  fpreal maxsubstep);
403 
404  /// Integrates the system for the given timestep using the
405  /// Euler method.
406  SIM_Result doEuler(SIM_Engine &engine,
407  SIM_ObjectArray &objects,
408  SIM_ObjectArray &newobjects,
409  SIM_ObjectArray &feedbacktoobjects,
410  const SIM_Time &starttime,
411  const SIM_Time &endtime,
412  const SIM_Time &timestep,
413  const SIM_Time &scaledtime);
414 
415  /// Integrates the system for the given timestep using the
416  /// Midpoint method.
417  SIM_Result doMidpoint(SIM_Engine &engine,
418  SIM_ObjectArray &objects,
419  SIM_ObjectArray &newobjects,
420  SIM_ObjectArray &feedbacktoobjects,
421  const SIM_Time &starttime,
422  const SIM_Time &endtime,
423  const SIM_Time &timestep,
424  const SIM_Time &scaledtime);
425 
426  /// Integrates the system for the given timestep using the
427  /// Runge-Kutta method.
428  SIM_Result doRungeKutta(SIM_Engine &engine,
429  SIM_ObjectArray &objects,
430  SIM_ObjectArray &newobjects,
431  SIM_ObjectArray &feedbacktoobjects,
432  const SIM_Time &starttime,
433  const SIM_Time &endtime,
434  const SIM_Time &timestep,
435  const SIM_Time &scaledtime);
436 
437  /// Integrate the system using the Runge-Kutta-Fehlberg method.
438  /// We also return a flag indicating if this failed, and a new
439  /// recommended substep value.
440  SIM_Result doRungeKuttaFehlberg(SIM_Engine &engine,
441  SIM_ObjectArray &objects,
442  SIM_ObjectArray &newobjects,
443  SIM_ObjectArray &feedbacktoobjects,
444  const SIM_Time &starttime,
445  const SIM_Time &endtime,
446  const SIM_Time &timestep,
447  const SIM_Time &scaledtime);
448 
449  /// Integrate the system using an embedded Runge-Kutta 2(1) method.
450  /// We also return a flag indicating if this failed, and a new
451  /// recommended substep value.
452  SIM_Result doRungeKutta21(SIM_Engine &engine,
453  SIM_ObjectArray &objects,
454  SIM_ObjectArray &newobjects,
455  SIM_ObjectArray &feedbacktoobjects,
456  const SIM_Time &starttime,
457  const SIM_Time &endtime,
458  const SIM_Time &timestep,
459  const SIM_Time &scaledtime);
460 
461  /// Integrate the system using an embedded Runge-Kutta 3(2) method.
462  /// We also return a flag indicating if this failed, and a new
463  /// recommended substep value.
464  SIM_Result doRungeKutta32(SIM_Engine &engine,
465  SIM_ObjectArray &objects,
466  SIM_ObjectArray &newobjects,
467  SIM_ObjectArray &feedbacktoobjects,
468  const SIM_Time &starttime,
469  const SIM_Time &endtime,
470  const SIM_Time &timestep,
471  const SIM_Time &scaledtime);
472 
473  /// Leapfrog integrator - a second order integrator requiring
474  /// only a single force evaluation per time step.
475  SIM_Result doLeapfrog(SIM_Engine &engine,
476  SIM_ObjectArray &objects,
477  SIM_ObjectArray &newobjects,
478  SIM_ObjectArray &feedbacktoobjects,
479  const SIM_Time &starttime,
480  const SIM_Time &endtime,
481  const SIM_Time &timestep,
482  const SIM_Time &scaledtime);
483 
484  // Helpfer functions for the Leapfrog method. The first updates
485  // particle positions according to the leapfrog method, while the
486  // second updates velocities.
487  void leapfrogAdvectObjs(SIM_ObjectArray &objects,
488  const SIM_Time &timestep);
489 
490  void leapfrogAdvect(SIM_Object *obj, GU_Detail *gdp,
491  const GA_PointGroup *ptgrp,
492  const SIM_Time &timestep);
493 
494  void leapfrogIntegrateObjs(SIM_ObjectArray &objects,
495  const SIM_Time &timestep);
496 
497  void leapfrogIntegrate(GU_Detail *gdp,
498  const GA_PointGroup *ptgrp,
499  const SIM_Time &timestep);
500 
501  /// Helper functions for the Leapfrog method. These apply
502  /// the position and velocity symplectic shear operators
503  /// to the system.
504  void leapfrogPositionUpdate(SIM_ObjectArray &objects,
505  const SIM_Time &timestep,
506  fpreal scale);
507 
508  void leapfrogPositionUpdate(SIM_Object *obj,
509  GU_Detail *gdp,
510  const GA_PointGroup *ptgrp,
511  const SIM_Time &timestep,
512  fpreal scale);
513 
514  void leapfrogVelocityUpdate(SIM_Engine &engine,
515  SIM_ObjectArray &objects,
516  SIM_ObjectArray &newobjects,
517  SIM_ObjectArray &feedbacktoobjects,
518  const SIM_Time &timestep,
519  fpreal scale);
520 
521 #if 0
522  /// FIXME
523  /// Integrate the system with the runge-kutta 4 method using
524  /// our new setup
525  ///
526  /// Note: solveRungeKuttaSteps now defaults to cplUpdate
527  /// instead of Taylor updates, this is now better
528  ///
529  void doNewRungeKutta(SIM_Engine &engine,
530  SIM_ObjectArray &objects,
531  SIM_ObjectArray &newobjects,
532  SIM_ObjectArray &feedbacktoobjects,
533  const SIM_Time &timestep);
534 #endif
535 
536  /// Helper function used to initialize our stored object info.
537  /// It is also used to intialize Runge Kutta data tables.
538  void initObjectInfoArray(SIM_ObjectArray &objects,
539  int numRungeKuttaSteps = 0);
540 
541  /// Applies all of our subsolvers to the current geometry. This
542  /// function exists so that this can easily be done multiple times
543  /// during a multi-step method.
544  void applyAllSolvers(SIM_Engine &engine,
545  SIM_ObjectArray &objects,
546  SIM_ObjectArray &newobjects,
547  SIM_ObjectArray &feedbacktoobjects,
548  const SIM_Time &timestep,
549  fpreal64 timescale,
550  bool incremental = false);
551 
552 #if 0
553  /// Applies the position solver to each object using the
554  /// given timestep.
555  void applyPositionSolver(SIM_Engine &engine,
556  SIM_ObjectArray &objects,
557  SIM_ObjectArray &newobjects,
558  SIM_ObjectArray &feedbacktoobjects,
559  const SIM_Time &timestep,
560  fpreal64 timescale);
561 #endif
562 
563  /// Clears the force attribute in all suplied data objects.
564  void clearAllForces(SIM_ObjectArray &objects);
565 
566  /// Integrates the current forces acting on the system for one
567  /// timestep using the standard Euler method.
568  /// currentscale: the amount by which to scale the deltas
569  /// when adding to the initial values for
570  /// the particle.
571  /// totalscale: the amount by which to scale the deltas
572  /// when adding them to the total velocity
573  /// change.
574  void applyCurrentDeltas(GU_Detail *gdp,
575  const GA_PointGroup *ptgrp,
576  const SIM_Time &timestep,
577  ObjectInfo &inf,
578  fpreal64 currentscale,
579  fpreal64 totalscale,
580  fpreal64 currentvscale,
581  fpreal64 totalvscale);
582 
583  /// Iterates through each object and calls the above function
584  /// on it.
585  void applyDeltasToObjs(SIM_ObjectArray &objects,
586  const SIM_Time &timestep,
587  fpreal64 currentscale,
588  fpreal64 totalscale,
589  fpreal64 currentvscale,
590  fpreal64 totalvscale);
591 
592  /// Applies changes to the current velocity of particles in the
593  /// geometry based on the contents of the force attribute. Also
594  /// applies the necessary changes to particle positions, depending
595  /// on the value of predict. If predict is false, then particle
596  /// positions are integrated based on the starting velocity (before
597  /// applying forces). Otherwise, they are integrated using the
598  /// post-force velocity.
599  void applyVelocityChanges(SIM_Object *obj, GU_Detail *gdp,
600  const GA_PointGroup *ptgrp,
601  const SIM_Time &timestep,
602  ObjectInfo &inf,
603  fpreal64 posscale = 1.0,
604  bool predict = true,
605  fpreal64 velscale = 1.0);
606 
607  void applyVelocityChangesToObjs(SIM_ObjectArray &objects,
608  const SIM_Time &timestep,
609  fpreal64 posscale = 1.0,
610  bool predict = true,
611  fpreal64 velscale = 1.0);
612 
613  THREADED_METHOD7(GAS_Integrator, gdp->getNumPoints() > 2048,
614  computeVelocityChanges,
615  GU_Detail *, gdp,
616  const GA_PointGroup *, ptgrp,
618  UT_Vector3Array &, velocities,
619  float, posscale,
620  bool, predict,
621  float, velscale)
622  void computeVelocityChangesPartial(
623  GU_Detail *gdp,
624  const GA_PointGroup *ptgrp,
625  SIM_Time timestep,
626  UT_Vector3Array &velocities,
627  float posscale,
628  bool predict,
629  float velscale,
630  const UT_JobInfo &info);
631 
632  /// Applies an incremental velocity change to all
633  /// objects, and adjusts the positions of particles
634  /// accordingly.
635  void applyVelocityIncrementsToObjs(SIM_ObjectArray &objects,
636  const SIM_Time &timestep);
637 
638  /// Advects the particles in the given geometry using the
639  /// given array of "velocity" vectors. We use an input
640  /// array here to allow for modifications to the velocity
641  /// (ie. to incorporate instantaneous changes due to acceleration
642  /// etc.). Using this method will allow us to support more
643  /// sophisticated advection techniques.
644  /// This method uses the AdvectionTYpe to determine how to act.
645  void advectParticles(SIM_Object *obj,
646  GU_Detail *gdp,
647  const GA_PointGroup *ptgrp,
648  const SIM_Time &timestep,
649  const UT_Vector3Array &velocities);
650 
651  /// Performs a simple pos += v * ts update
652  THREADED_METHOD4(GAS_Integrator, gdp->getNumPoints() > 2048,
653  advectParticlesEuler,
654  GU_Detail *, gdp,
655  const GA_PointGroup *, ptgrp,
657  const UT_Vector3Array &, velocities)
658  void advectParticlesEulerPartial(GU_Detail *gdp,
659  const GA_PointGroup *ptgrp,
660  SIM_Time timestep,
661  const UT_Vector3Array &velocities,
662  const UT_JobInfo &info);
663 
664  /// Performs a euler update on the attributes
665  /// Simplistic Euler update with no temp arrays.
666  void directIntegratePosObjs(SIM_ObjectArray &objects,
667  const SIM_Time &timestep);
668  THREADED_METHOD3(GAS_Integrator, gdp->getNumPoints() > 2048,
669  directIntegratePos,
670  GU_Detail *, gdp,
671  const GA_PointGroup *, ptgrp,
673  void directIntegratePosPartial(GU_Detail *gdp,
674  const GA_PointGroup *ptgrp,
675  SIM_Time timestep,
676  const UT_JobInfo &info);
677 
678  THREADED_METHOD3(GAS_Integrator, gdp->getNumPoints() > 2048,
679  directIntegrateOrient,
680  GU_Detail *, gdp,
681  const GA_PointGroup *, ptgrp,
683  void directIntegrateOrientPartial(GU_Detail *gdp,
684  const GA_PointGroup *ptgrp,
685  SIM_Time timestep,
686  const UT_JobInfo &info);
687 
688  /// Computes the xsph velocity given the specified velocity array
689  THREADED_METHOD6(GAS_Integrator, gdp->getNumPoints() > 2048,
690  computeXSPHVel,
691  GU_Detail *, gdp,
692  const GA_PointGroup *, ptgrp,
693  const GAS_SPH &, sph,
694  const UT_Vector3Array &, refvel,
695  UT_Vector3Array &, xsphvel,
696  float, xsphConstant)
697  void computeXSPHVelPartial(GU_Detail *gdp,
698  const GA_PointGroup *ptgrp,
699  const GAS_SPH &sph,
700  const UT_Vector3Array &refvel,
701  UT_Vector3Array &xsphvel,
702  float xsphConstant,
703  const UT_JobInfo &info);
704 
705  /// Resets positions and velocities back to their current values.
706  void resetObjectStates(SIM_ObjectArray &objects);
707 
708  /// Sets the final value of velocity to the total velocity
709  /// accumulated for each point, and sets each point's position
710  /// back to its initial value.
711  void setFinalValues(SIM_ObjectArray &objects);
712 
713  /// Responds to impact information attached to a particle by
714  /// bouncing the particle from the impact point.
715  /// Returns true if anything was changed.
716  bool bounceParticleFromImpact(GU_Detail *gdp,
717  const SIM_Time &starttime,
718  const SIM_Time &endtime,
719  const SIM_Time &timestep,
720  const UT_DMatrix4 &xform,
721  const UT_DMatrix4 &ixform,
722  const SIM_Impacts *impacts,
723  int impactidx, GA_Index ptnum);
724 
725  /// Responds to all impacts for a given object.
726  void respondToCollisions(SIM_Engine &engine,
727  SIM_Object &object,
728  const SIM_ObjectArray &affectors,
729  const UT_DMatrix4 &xform,
730  const UT_DMatrix4 &ixform,
731  const UT_DMatrix4 &posxform,
732  const SIM_Time &starttime,
733  const SIM_Time &endtime,
734  const SIM_Time &timestep,
735  const char *impactdataname,
736  bool clearimpacts = false);
737 
738  /// Applies impacts to the particles stored in object.
739  void handleObjectImpacts(SIM_Engine &engine,
740  SIM_Object &object,
741  const SIM_Time &starttime,
742  const SIM_Time &endtime,
743  const SIM_Time &timestep,
744  bool clearimpacts = false);
745 
746  /// Applies impacts to all objects.
747  void handleAllImpacts(SIM_Engine &engine,
748  SIM_ObjectArray &objects,
749  const SIM_Time &starttime,
750  const SIM_Time &endtime,
751  const SIM_Time &timestep,
752  bool clearimpacts = false);
753 
754  /// Do any post processing after a timestep
755  void doPostSolveUpdates( SIM_Engine &engine,
756  SIM_ObjectArray &objects,
757  SIM_ObjectArray &newobjects,
758  SIM_ObjectArray &feedbacktoobjects,
759  const SIM_Time &timestep );
760 
761  static const SIM_DopDescription *getDopDescription();
762 
763  /// Gets the solver to use for getting object properties.
764  const SIM_Solver *getCurrentSolver() const;
765 
766  /// We store an array of object info objects to keep track of any
767  /// intermediate object states during multi-step integration
768  /// methods.
769  UT_Array<ObjectInfo> myObjectInfoArray;
770 
771  /// We store a persistent array of runge-kutta data tables, each
772  /// of which stores intermediate RK integration steps for a
773  /// different object. We store this permanently so that we don't
774  /// have to repeatedly allocate the (undoubedly large amount of)
775  /// memory required for this data structure on each substep.
776  UT_Array<RungeKuttaDataTable> myRungeKuttaData;
777 
778  /// Store a pointer to an engine so that we can add errors
779  /// as necessary.
780  SIM_Engine *myEngine;
781 
782  // Store matrices of values for each type of Runge-Kutta integration
783  // we want to do
784  static const fpreal64 RK_45_COEFS[][GAS_RK_45_NUMSTEPS-1];
785  static const fpreal64 RK_45_TIME_COEFS[GAS_RK_45_NUMSTEPS-1];
786  static const fpreal64 RK_45_TOTAL_COEFS4[GAS_RK_45_NUMSTEPS];
787  static const fpreal64 RK_45_TOTAL_COEFS5[GAS_RK_45_NUMSTEPS];
788 
789  static const fpreal64 RK_4_COEFS[][GAS_RK_4_NUMSTEPS-1];
790  static const fpreal64 RK_4_TIME_COEFS[GAS_RK_4_NUMSTEPS-1];
791  static const fpreal64 RK_4_TOTAL_COEFS[GAS_RK_4_NUMSTEPS];
792 
793  static const fpreal64 RK_21_COEFS[][GAS_RK_21_NUMSTEPS-1];
794  static const fpreal64 RK_21_TIME_COEFS[GAS_RK_21_NUMSTEPS-1];
795  static const fpreal64 RK_21_TOTAL_COEFS1[GAS_RK_21_NUMSTEPS];
796  static const fpreal64 RK_21_TOTAL_COEFS2[GAS_RK_21_NUMSTEPS];
797  static const fpreal64 RK_21_TOLERANCE_COEF;
798 
799  static const fpreal64 RK_32_COEFS[][GAS_RK_32_NUMSTEPS-1];
800  static const fpreal64 RK_32_TIME_COEFS[GAS_RK_32_NUMSTEPS-1];
801  static const fpreal64 RK_32_TOTAL_COEFS2[GAS_RK_32_NUMSTEPS];
802  static const fpreal64 RK_32_TOTAL_COEFS3[GAS_RK_32_NUMSTEPS];
803 
804  static const fpreal64 EULER_TOTAL_COEFS[GAS_EULER_NUMSTEPS];
805 
806  static const fpreal64 MIDPOINT_COEFS[][GAS_MIDPOINT_NUMSTEPS-1];
807  static const fpreal64 MIDPOINT_TIME_COEFS[GAS_MIDPOINT_NUMSTEPS-1];
808  static const fpreal64 MIDPOINT_TOTAL_COEFS[GAS_MIDPOINT_NUMSTEPS];
809 
810  SIM_Solver *myCurrentSolver;
811 
815  "GAS Integrator",
816  getDopDescription());
817 };
818 
819 #endif
820 
#define DECLARE_STANDARD_GETCASTTOTYPE()
Definition: SIM_DataUtils.h:45
#define THREADED_METHOD7(CLASSNAME, DOMULTI, METHOD, PARMTYPE1, PARMNAME1, PARMTYPE2, PARMNAME2, PARMTYPE3, PARMNAME3, PARMTYPE4, PARMNAME4, PARMTYPE5, PARMNAME5, PARMTYPE6, PARMNAME6, PARMTYPE7, PARMNAME7)
#define SIM_NAME_PRIMARYSOLVER
Definition: SIM_Names.h:169
const GA_PointGroup * ptgrp
#define GAS_NAME_GEOMETRY
Definition: GAS_Utils.h:28
#define GAS_NAME_USETIMESTEP
Definition: GAS_Utils.h:37
#define GAS_API
Definition: GAS_API.h:10
3D Vector class.
png_uint_32 i
Definition: png.h:2877
virtual bool solveGasSubclass(SIM_Engine &engine, SIM_Object *obj, SIM_Time time, SIM_Time timestep)
exint GA_Size
Defines the bit width for index and offset types in GA.
Definition: GA_Types.h:211
#define GAS_RK_32_NUMSTEPS
#define THREADED_METHOD3(CLASSNAME, DOMULTI, METHOD, PARMTYPE1, PARMNAME1, PARMTYPE2, PARMNAME2, PARMTYPE3, PARMNAME3)
#define DECLARE_DATAFACTORY(DataClass, SuperClass, Description, DopParms)
Definition: SIM_DataUtils.h:58
GA_API const UT_StringHolder scale
#define GET_DATA_FUNC_I(DataName, FuncName)
#define GAS_MIDPOINT_NUMSTEPS
Holds pointers to a number of SIM_Object objects.
double fpreal64
Definition: SYS_Types.h:185
virtual void getImpulseMassMatrixSubclass(const SIM_Object &object, const UT_Vector3 &impulseworldpos, UT_DMatrix3 &immatrix) const
#define GETSET_DATA_FUNCS_I(DataName, FuncName)
virtual SIM_PointImpulseMassMatrixResolver * getPointImpulseMassMatrixResolverSubclass(const SIM_Object &object) const
Builds a resolver for evaluating mass matrices swiftly.
#define THREADED_METHOD5(CLASSNAME, DOMULTI, METHOD, PARMTYPE1, PARMNAME1, PARMTYPE2, PARMNAME2, PARMTYPE3, PARMNAME3, PARMTYPE4, PARMNAME4, PARMTYPE5, PARMNAME5)
GLboolean * data
Definition: glcorearb.h:130
#define THREADED_METHOD4(CLASSNAME, DOMULTI, METHOD, PARMTYPE1, PARMNAME1, PARMTYPE2, PARMNAME2, PARMTYPE3, PARMNAME3, PARMTYPE4, PARMNAME4)
#define GAS_RK_45_NUMSTEPS
GA_Size GA_Index
Define the strictness of GA_Offset/GA_Index.
Definition: GA_Types.h:611
virtual SIM_Result solveObjectsSubclass(SIM_Engine &engine, SIM_ObjectArray &objects, SIM_ObjectArray &newobjects, SIM_ObjectArray &feedbacktoobjects, const SIM_Time &timestep)
Merely calls solve on each object.
#define GAS_RK_4_NUMSTEPS
double fpreal
Definition: SYS_Types.h:263
typedef int
Definition: png.h:1175
#define GET_DATA_FUNC_B(DataName, FuncName)
#define GET_DATA_FUNC_F(DataName, FuncName)
#define THREADED_METHOD6(CLASSNAME, DOMULTI, METHOD, PARMTYPE1, PARMNAME1, PARMTYPE2, PARMNAME2, PARMTYPE3, PARMNAME3, PARMTYPE4, PARMNAME4, PARMTYPE5, PARMNAME5, PARMTYPE6, PARMNAME6)
#define GAS_NAME_TIMESCALE
Definition: GAS_Utils.h:38
bool shouldMultiThread()
#define GAS_EULER_NUMSTEPS
virtual void getPointImpulseMassMatrixSubclass(const SIM_Object &object, int ptnum, UT_DMatrix3 &immatrix) const
#define GET_DATA_FUNC_S(DataName, FuncName)
#define GAS_RK_21_NUMSTEPS
const GA_PointGroup SIM_Time timestep
This implements a SIM_Geometry that copies the source geometry.