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