initial commit of the fluid simulator.
[blender.git] / intern / elbeem / intern / simulation_object.cpp
1 /******************************************************************************
2  *
3  * El'Beem - Free Surface Fluid Simulation with the Lattice Boltzmann Method
4  * Copyright 2003,2004 Nils Thuerey
5  *
6  * Basic interface for all simulation modules
7  *
8  *****************************************************************************/
9
10 #include "simulation_object.h"
11 #include "ntl_bsptree.h"
12 #include "ntl_scene.h"
13 #include "factory_lbm.h"
14 #include "lbmfunctions.h"
15
16 #ifdef _WIN32
17 #else
18 #include <sys/time.h>
19 #endif
20
21
22
23
24 /******************************************************************************
25  * Constructor
26  *****************************************************************************/
27 SimulationObject::SimulationObject() :
28         ntlGeometryShader(),
29         mGeoStart(-100.0), mGeoEnd(100.0),
30         mGeoInitId(-1), mpGiTree(NULL), mpGiObjects(NULL),
31         mpGlob(NULL),
32         mPanic( false ),
33         mDebugType( 1 /* =FLUIDDISPNothing*/ ),
34         mSolverType("-"), mStepsPerFrame( 10 ),
35         mpLbm(NULL),
36         mpParam( NULL ),
37         mShowSurface(true), mShowParticles(false),
38         mSelectedCid( NULL ),
39
40         stnOld("opt"),
41         stnFsgr("fsgr")
42 {
43         mpParam = new Parametrizer();
44
45         for(int i=0; i<MAX_DEBDISPSET; i++) {
46                 mDebDispSet[i].type  = (i);
47                 mDebDispSet[i].on    = false;
48                 mDebDispSet[i].scale = 1.0;
49         }
50
51         // reset time
52         mTime                                           = 0.0;
53         mDisplayTime                    = 0.0;
54 }
55
56
57 /******************************************************************************
58  * Destructor
59  *****************************************************************************/
60 SimulationObject::~SimulationObject()
61 {
62         if(mpGiTree != NULL) delete mpGiTree;
63         delete mpLbm;
64   delete mpParam;
65         debMsgStd("SimulationObject",DM_MSG,"El'Beem Done!\n",10);
66 }
67
68
69
70 /*****************************************************************************/
71 /*! init tree for certain geometry init */
72 /*****************************************************************************/
73 void SimulationObject::initGeoTree(int id) {
74         if(mpGlob == NULL) { errorOut("SimulationObject::initGeoTree error: Requires globals!"); exit(1); }
75         mGeoInitId = id;
76         ntlScene *scene = mpGlob->getScene();
77         mpGiObjects = scene->getObjects();
78
79         if(mpGiTree != NULL) delete mpGiTree;
80         char treeFlag = (1<<(mGeoInitId+4));
81         mpGiTree = new ntlTree( 20, 4, // warning - fixed values for depth & maxtriangles here...
82                                                                                                 scene, treeFlag );
83 }
84
85 /*****************************************************************************/
86 /*! destroy tree etc. when geometry init done */
87 /*****************************************************************************/
88 void SimulationObject::freeGeoTree() {
89         if(mpGiTree != NULL) delete mpGiTree;
90 }
91
92
93
94
95 /******************************************************************************
96  * simluation interface: initialize simulation using the given configuration file 
97  *****************************************************************************/
98 int SimulationObject::initializeLbmSimulation(ntlRenderGlobals *glob)
99 {
100         mpGlob = glob;
101         if(!getVisible()) {
102                 mpAttrs->setAllUsed();
103                 return 0;
104         }
105
106         //mDimension is deprecated
107         mSolverType = mpAttrs->readString("solver", mSolverType, "SimulationObject","mSolverType", false ); 
108         if(mSolverType == stnFsgr) {
109                 mpLbm = createSolverLbmFsgr(); 
110         } else if(mSolverType == stnOld) {
111 #ifdef LBM_INCLUDE_TESTSOLVERS
112                 // old solver for gfx demo
113                 mpLbm = createSolverOld();
114 #endif // LBM_TESTSOLVER
115         } else {
116                 errorOut("SimulationObject::initializeLbmSimulation : Invalid solver type - note that mDimension is deprecated, use the 'solver' keyword instead");
117                 exit(1);
118         }
119
120
121   /* check lbm pointer */
122         if(mpLbm == NULL) {
123                 errorOut("SimulationObject::initializeLbmSimulation : Unable to init dim"<<mSolverType<<" LBM solver! ");
124                 exit(1);
125         }
126         debugOut("SimulationObject::initialized "<< mpLbm->getIdString() <<" LBM solver! ", 2);
127
128         // for non-param simulations
129         mpLbm->setParametrizer( mpParam );
130         mpParam->setAttrList( getAttributeList() );
131         mpParam->setSize( mpLbm->getSizeX(), mpLbm->getSizeY(), mpLbm->getSizeZ() );
132         mpParam->parseAttrList();
133
134         mpLbm->setAttrList( getAttributeList() );
135         mpLbm->parseAttrList();
136         mParts.parseAttrList( getAttributeList() );
137         mParts.setName( getName() + "_part" );
138         mParts.initialize( glob );
139         
140         // init material settings
141         string matMc("default");
142         matMc = mpAttrs->readString("material_surf", matMc, "SimulationObject","matMc", false );
143         mShowSurface   = mpAttrs->readInt("showsurface", mShowSurface, "SimulationObject","mShowSurface", false ); 
144         mShowParticles = mpAttrs->readInt("showparticles", mShowParticles, "SimulationObject","mShowParticles", false ); 
145
146         checkBoundingBox( mGeoStart, mGeoEnd, "SimulationObject::initializeSimulation" );
147         mpLbm->setGeoStart( mGeoStart );
148         mpLbm->setGeoEnd( mGeoEnd );
149         mpLbm->setRenderGlobals( mpGlob );
150         mpLbm->setName( getName() + "_lbm" );
151         mpLbm->initialize( NULL, mpGlob->getScene()->getObjects() );
152
153         // print cell type stats
154         const int jmax = sizeof(CellFlagType)*8;
155         int totalCells = 0;
156         int flagCount[jmax];
157         for(int j=0; j<jmax ; j++) flagCount[j] = 0;
158         int diffInits = 0;
159         LbmSolverInterface::CellIdentifier cid = mpLbm->getFirstCell();
160         for(; mpLbm->noEndCell( cid );
161               mpLbm->advanceCell( cid ) ) {
162                 int flag = mpLbm->getCellFlag(cid,0);
163                 int flag2 = mpLbm->getCellFlag(cid,1);
164                 if(flag != flag2) {
165                         diffInits++;
166                 }
167                 for(int j=0; j<jmax ; j++) {
168                         if( flag&(1<<j) ) flagCount[j]++;
169                 }
170                 totalCells++;
171         }
172         mpLbm->deleteCellIterator( &cid );
173
174 #if ELBEEM_BLENDER!=1
175         char charNl = '\n';
176         debugOutNnl("SimulationObject::initializeLbmSimulation celltype stats: " <<charNl, 5);
177         debugOutNnl("no. of cells = "<<totalCells<<", "<<charNl ,5);
178         for(int j=0; j<jmax ; j++) {
179                 std::ostringstream out;
180                 if(flagCount[j]>0) {
181                         out<<"\t" << flagCount[j] <<" x "<< convertCellFlagType2String( (CellFlagType)(1<<j) ) <<", " << charNl;
182                         debugOutNnl(out.str(), 5);
183                 }
184         }
185         // compute dist. of empty/bnd - fluid - if
186         // cfEmpty   = (1<<0), cfBnd  = (1<< 2), cfFluid   = (1<<10), cfInter   = (1<<11),
187         {
188                 std::ostringstream out;
189                 out.precision(2); out.width(4);
190                 int totNum = flagCount[0]+flagCount[2]+flagCount[10]+flagCount[11];
191                 double ebFrac = (double)(flagCount[0]+flagCount[2]) / totNum;
192                 double flFrac = (double)(flagCount[10]) / totNum;
193                 double ifFrac = (double)(flagCount[11]) / totNum;
194                 double eifFrac = (double)(flagCount[11]+flagCount[26]+flagCount[27]) / totNum;
195                 //???
196                 out<<"\tFractions: [empty/bnd - fluid - interface - ext. if]  =  [" << ebFrac<<" - " << flFrac<<" - " << ifFrac<<" - " << eifFrac <<"] "<< charNl;
197
198                 if(diffInits > 0) {
199                         debugOut("SimulationObject::initializeLbmSimulation celltype Warning: Diffinits="<<diffInits<<" !!!!!!!!!" , 5);
200                 }
201                 debugOutNnl(out.str(), 5);
202         }
203 #endif // ELBEEM_BLENDER==1
204
205         mpLbm->initParticles( &mParts );
206
207         // this has to be inited here - before, the values might be unknown
208         ntlGeometryObject *surf = mpLbm->getSurfaceGeoObj();
209         if(surf) {
210                 surf->setName( "final" ); // final surface mesh 
211                 // warning - this might cause overwriting effects for multiple sims and geom dump...
212                 surf->setCastShadows( true );
213                 surf->setReceiveShadows( false );
214                 surf->searchMaterial( glob->getMaterials() );
215                 if(mShowSurface) mObjects.push_back( surf );
216         }
217         
218         mParts.setStart( mGeoStart );
219         mParts.setEnd( mGeoEnd );
220         mParts.setCastShadows( false );
221         mParts.setReceiveShadows( false );
222         mParts.searchMaterial( glob->getMaterials() );
223         if(mShowParticles) mObjects.push_back( &mParts );
224
225         // add objects to display for debugging (e.g. levelset particles)
226         vector<ntlGeometryObject *> debugObjs = mpLbm->getDebugObjects();
227         for(size_t i=0;i<debugObjs.size(); i++) {
228                 debugObjs[i]->setCastShadows( false );
229                 debugObjs[i]->setReceiveShadows( false );
230                 debugObjs[i]->searchMaterial( glob->getMaterials() );
231                 mObjects.push_back( debugObjs[i] );
232         }
233         return 0;
234 }
235
236
237 /******************************************************************************
238  * simluation interface: advance simulation another step (whatever delta time that might be) 
239  *****************************************************************************/
240 void SimulationObject::step( void )
241 {
242         mpLbm->step();
243
244         mParts.savePreviousPositions();
245         mpLbm->advanceParticles( &mParts );
246         mTime += mpParam->getStepTime();
247         if(mpLbm->getPanic()) mPanic = true;
248
249   //debMsgStd("SimulationObject::step",DM_MSG," Sim '"<<mName<<"' stepped to "<<mTime<<" (stept="<<(mpParam->getStepTime())<<", framet="<<getFrameTime()<<") ", 10);
250 }
251 /*! prepare visualization of simulation for e.g. raytracing */
252 void SimulationObject::prepareVisualization( void ) {
253         if(mPanic) return;
254         mpLbm->prepareVisualization();
255 }
256
257
258 /******************************************************************************/
259 /* get current start simulation time */
260 double SimulationObject::getStartTime( void ) {
261         //return mpParam->calculateAniStart();
262         return mpParam->getAniStart();
263 }
264 /* get time for a single animation frame */
265 double SimulationObject::getFrameTime( void ) {
266         return mpParam->getAniFrameTime();
267 }
268 /* get time for a single time step  */
269 double SimulationObject::getStepTime( void ) {
270         return mpParam->getStepTime();
271 }
272
273
274 /******************************************************************************
275  * return a pointer to the geometry object of this simulation 
276  *****************************************************************************/
277 //ntlGeometryObject *SimulationObject::getGeometry() { return mpMC; }
278 vector<ntlGeometryObject *>::iterator 
279 SimulationObject::getObjectsBegin()
280 {
281         return mObjects.begin();
282 }
283 vector<ntlGeometryObject *>::iterator 
284 SimulationObject::getObjectsEnd()
285 {
286         return mObjects.end();
287 }
288
289
290
291
292
293 /******************************************************************************
294  * GUI - display debug info 
295  *****************************************************************************/
296
297 void SimulationObject::drawDebugDisplay() {
298 #ifndef NOGUI
299         //debugOut(" SD: "<<mDebugType<<" v"<<getVisible()<<" don"<< (mDebDispSet[mDebugType].on) , 10);
300         if(!getVisible()) return;
301
302         if( mDebugType > (MAX_DEBDISPSET-1) ){
303                 errorOut("SimulationObject::drawDebugDisplay : Invalid debug type!");
304                 exit(1);
305         }
306
307         mDebDispSet[ mDebugType ].on = true;
308         //errorOut( mDebugType <<"//"<< mDebDispSet[mDebugType].type );
309         mpLbm->debugDisplay( &mDebDispSet[ mDebugType ] );
310
311         ::lbmMarkedCellDisplay<>( mpLbm );
312 #endif
313 }
314
315 /* GUI - display interactive info  */
316 void SimulationObject::drawInteractiveDisplay()
317 {
318 #ifndef NOGUI
319         if(!getVisible()) return;
320         if(mSelectedCid) {
321                 // in debugDisplayNode if dispset is on is ignored...
322                 ::debugDisplayNode<>( &mDebDispSet[ FLUIDDISPGrid ], mpLbm, mSelectedCid );
323         }
324 #endif
325 }
326
327
328 /*******************************************************************************/
329 // GUI - handle mouse movement for selection 
330 /*******************************************************************************/
331 void SimulationObject::setMousePos(int x,int y, ntlVec3Gfx org, ntlVec3Gfx dir)
332 {
333         normalize( dir );
334         // assume 2D sim is in XY plane...
335         
336         double zplane = (mGeoEnd[2]-mGeoStart[2])*0.5;
337         double zt = (zplane-org[2]) / dir[2];
338         ntlVec3Gfx pos(
339                         org[0]+ dir[0] * zt,
340                         org[1]+ dir[1] * zt, 0.0);
341
342         mSelectedCid = mpLbm->getCellAt( pos );
343         //errMsg("SMP ", mName<< x<<" "<<y<<" - "<<dir );
344         x = y = 0; // remove warning
345 }
346                         
347
348 void SimulationObject::setMouseClick()
349 {
350         if(mSelectedCid) {
351                 ::debugPrintNodeInfo<>( mpLbm, mSelectedCid, mpLbm->getNodeInfoString() );
352         }
353 }
354
355