Several minor fixes:
authorNils Thuerey <nils@thuerey.de>
Wed, 29 Mar 2006 07:35:54 +0000 (07:35 +0000)
committerNils Thuerey <nils@thuerey.de>
Wed, 29 Mar 2006 07:35:54 +0000 (07:35 +0000)
- Added part of Austin's msvc8 fixes (vector::erase function
  was "misused"), hopefully compiles better now.
- Ctrl-b now also bakes a selected fluidsim domain
  similar to the softbodies.
- Added surface smoothing option for domains: default is
  1, higher values result in a smoother surface (and probably
  slightly higher comupation times), while 0 means the surface
  is not modified at all.
- Added BLENDER_ELBEEMBOBJABORT environment variable in readBobj,
  if >0 quits blender when a not yet existing fluidsim
  frame should be loaded. Useful for rendering simulations
  as far as possible from the command line.
- Surface normals pointer is now set to NULL in readfile.c
- Fixed win32 error string handling, now uses a function
  to return the string from the solver.
- Fixed fluidsim particle halo scaling problem.
- Solver update

30 files changed:
intern/elbeem/extern/LBM_fluidsim.h
intern/elbeem/extern/elbeem.h
intern/elbeem/intern/elbeem.cpp
intern/elbeem/intern/elbeem.h
intern/elbeem/intern/isosurface.cpp
intern/elbeem/intern/isosurface.h
intern/elbeem/intern/ntl_blenderdumper.cpp
intern/elbeem/intern/ntl_geometryobject.cpp
intern/elbeem/intern/ntl_geometryobject.h
intern/elbeem/intern/ntl_geometryshader.h
intern/elbeem/intern/ntl_world.h
intern/elbeem/intern/particletracer.cpp
intern/elbeem/intern/particletracer.h
intern/elbeem/intern/simulation_object.cpp
intern/elbeem/intern/simulation_object.h
intern/elbeem/intern/solver_class.h
intern/elbeem/intern/solver_init.cpp
intern/elbeem/intern/solver_interface.cpp
intern/elbeem/intern/solver_interface.h
intern/elbeem/intern/solver_main.cpp
intern/elbeem/intern/solver_relax.h
intern/elbeem/intern/solver_util.cpp
source/blender/blenkernel/intern/DerivedMesh.c
source/blender/blenkernel/intern/effect.c
source/blender/blenloader/intern/readfile.c
source/blender/makesdna/DNA_object_fluidsim.h
source/blender/render/intern/source/convertblender.c
source/blender/src/buttons_object.c
source/blender/src/fluidsim.c
source/blender/src/space.c

index b879dd4..c93d129 100644 (file)
@@ -68,16 +68,6 @@ int performElbeemSimulation(char *cfgfilename);
 void fluidsimGetAxisAlignedBB(struct Mesh *mesh, float obmat[][4],
                 /*RET*/ float start[3], /*RET*/ float size[3], /*RET*/ struct Mesh **bbmesh );
 
-// implemented in intern/elbeem/utilities.cpp
-/* set elbeem debug output level (0=off to 10=full on) */
-void elbeemSetDebugLevel(int level);
-/* elbeem debug output function */
-void elbeemDebugOut(char *msg);
-
-/* estimate how much memory a given setup will require */
-double elbeemEstimateMemreq(int res, 
-               float sx, float sy, float sz,
-               int refine, char *retstr);
 
 #endif
 
index 83a9cb9..f0d154c 100644 (file)
@@ -73,6 +73,10 @@ typedef struct elbeemSimulationSettings {
 
        /* global transformation to apply to fluidsim mesh */
        float surfaceTrafo[4*4];
+
+       /* development variables, testing for upcoming releases...*/
+       float farFieldSize;
+
 } elbeemSimulationSettings;
 
 
@@ -122,8 +126,12 @@ extern "C"  {
 // reset elbeemSimulationSettings struct with defaults
 void elbeemResetSettings(struct elbeemSimulationSettings*);
  
-// start fluidsim init
+// start fluidsim init (returns !=0 upon failure)
 int elbeemInit(struct elbeemSimulationSettings*);
+// get failure message during simulation or init
+// if an error occured (the string is copied into buffer,
+// max. length = 256 chars )
+void elbeemGetErrorString(char *buffer);
 
 // reset elbeemMesh struct with zeroes
 void elbeemResetMesh(struct elbeemMesh*);
@@ -135,17 +143,36 @@ int elbeemAddMesh(struct elbeemMesh*);
 int elbeemSimulate(void);
 
 
-// helper function - simplify animation channels
+// helper functions 
+
+// simplify animation channels
 // returns if the channel and its size changed
 int elbeemSimplifyChannelFloat(float *channel, int *size);
 int elbeemSimplifyChannelVec3(float *channel, int *size);
 
+// helper functions implemented in utilities.cpp
+
+/* set elbeem debug output level (0=off to 10=full on) */
+void elbeemSetDebugLevel(int level);
+/* elbeem debug output function, prints if debug level >0 */
+void elbeemDebugOut(char *msg);
+
+/* estimate how much memory a given setup will require */
+double elbeemEstimateMemreq(int res,
+    float sx, float sy, float sz,
+    int refine, char *retstr);
+
+
+
 #ifdef __cplusplus
 }
 #endif // __cplusplus
 
+
+
 /******************************************************************************/
-// internal defines, do not use for setting up simulation
+// internal defines, do not use for initializing elbeemMesh
+// structs, for these use OB_xxx defines above
 
 /*! fluid geometry init types */
 #define FGI_FLAGSTART   16
index 8e96879..daee5f8 100644 (file)
@@ -67,6 +67,8 @@ void elbeemResetSettings(elbeemSimulationSettings *set) {
        set->generateVertexVectors = 0;
        set->surfaceSmoothing = 1.;
 
+       set->farFieldSize = 0.;
+
        // init identity
        for(int i=0; i<16; i++) set->surfaceTrafo[i] = 0.0;
        for(int i=0; i<4; i++) set->surfaceTrafo[i*4+i] = 1.0;
@@ -87,6 +89,13 @@ int elbeemInit(elbeemSimulationSettings *settings) {
        return 0;
 }
 
+// error message access
+extern "C" 
+void elbeemGetErrorString(char *buffer) {
+       if(!buffer) return;
+       strncpy(buffer,gElbeemErrorString,256);
+}
+
 // reset elbeemMesh struct with zeroes
 extern "C" 
 void elbeemResetMesh(elbeemMesh *mesh) {
index 83a9cb9..f0d154c 100644 (file)
@@ -73,6 +73,10 @@ typedef struct elbeemSimulationSettings {
 
        /* global transformation to apply to fluidsim mesh */
        float surfaceTrafo[4*4];
+
+       /* development variables, testing for upcoming releases...*/
+       float farFieldSize;
+
 } elbeemSimulationSettings;
 
 
@@ -122,8 +126,12 @@ extern "C"  {
 // reset elbeemSimulationSettings struct with defaults
 void elbeemResetSettings(struct elbeemSimulationSettings*);
  
-// start fluidsim init
+// start fluidsim init (returns !=0 upon failure)
 int elbeemInit(struct elbeemSimulationSettings*);
+// get failure message during simulation or init
+// if an error occured (the string is copied into buffer,
+// max. length = 256 chars )
+void elbeemGetErrorString(char *buffer);
 
 // reset elbeemMesh struct with zeroes
 void elbeemResetMesh(struct elbeemMesh*);
@@ -135,17 +143,36 @@ int elbeemAddMesh(struct elbeemMesh*);
 int elbeemSimulate(void);
 
 
-// helper function - simplify animation channels
+// helper functions 
+
+// simplify animation channels
 // returns if the channel and its size changed
 int elbeemSimplifyChannelFloat(float *channel, int *size);
 int elbeemSimplifyChannelVec3(float *channel, int *size);
 
+// helper functions implemented in utilities.cpp
+
+/* set elbeem debug output level (0=off to 10=full on) */
+void elbeemSetDebugLevel(int level);
+/* elbeem debug output function, prints if debug level >0 */
+void elbeemDebugOut(char *msg);
+
+/* estimate how much memory a given setup will require */
+double elbeemEstimateMemreq(int res,
+    float sx, float sy, float sz,
+    int refine, char *retstr);
+
+
+
 #ifdef __cplusplus
 }
 #endif // __cplusplus
 
+
+
 /******************************************************************************/
-// internal defines, do not use for setting up simulation
+// internal defines, do not use for initializing elbeemMesh
+// structs, for these use OB_xxx defines above
 
 /*! fluid geometry init types */
 #define FGI_FLAGSTART   16
index ed0fab2..07bc6b7 100644 (file)
@@ -37,7 +37,7 @@ IsoSurface::IsoSurface(double iso) :
   mInitDone(false),
        mSmoothSurface(0.0), mSmoothNormals(0.0),
        mAcrossEdge(), mAdjacentFaces(),
-       mCutoff(-1), // off by default
+       mCutoff(-1), mCutArray(NULL),// off by default
        mFlagCnt(1),
        mSCrad1(0.), mSCrad2(0.), mSCcenter(0.)
 {
@@ -152,6 +152,7 @@ void IsoSurface::triangulate( void )
        const int cubieOffsetZ[8] = {
                0,0,0,0,  1,1,1,1 };
 
+       const int coAdd=2;
   // let the cubes march 
        pz = mStart[2]-gsz*0.5;
        for(int k=1;k<(mSizez-2);k++) {
@@ -259,12 +260,17 @@ void IsoSurface::triangulate( void )
 
                                }
 
-                               const int coAdd=2;
-                               if(i<coAdd+mCutoff) continue;
-                               if(j<coAdd+mCutoff) continue;
-                               if((mCutoff>0) && (k<coAdd)) continue;
-                               if(i>mSizex-2-coAdd-mCutoff) continue;
-                               if(j>mSizey-2-coAdd-mCutoff) continue;
+                               if( (i<coAdd+mCutoff) ||
+                                   (j<coAdd+mCutoff) ||
+                                   ((mCutoff>0) && (k<coAdd)) ||// bottom layer
+                                   (i>mSizex-2-coAdd-mCutoff) ||
+                                   (j>mSizey-2-coAdd-mCutoff) ) {
+                                       if(mCutArray) {
+                                               if(k < mCutArray[j*this->mSizex+i]) continue;
+                                       } else {
+                                               continue;
+                                       }
+                               }
 
                                // Create the triangles... 
                                for(int e=0; mcTriTable[cubeIndex][e]!=-1; e+=3) {
index 90a4e11..5e5115e 100644 (file)
@@ -95,6 +95,8 @@ class IsoSurface :
 
                //! cutoff border area
                int mCutoff;
+               //! cutoff heigh values
+               int *mCutArray;
                
                //! trimesh vars
                vector<int> flags;
@@ -160,6 +162,8 @@ class IsoSurface :
                }
                //! set cut off border
                inline void setCutoff(int set) { mCutoff = set; };
+               //! set cut off border
+               inline void setCutArray(int *set) { mCutArray = set; };
 
                //! OpenGL viz "interface"
                unsigned int getIsoVertexCount() {
index fce5a08..ecd1967 100644 (file)
@@ -91,7 +91,7 @@ int ntlBlenderDumper::renderScene( void )
   vector<ntlTriangle> Triangles;
   vector<ntlVec3Gfx>  Vertices;
   vector<ntlVec3Gfx>  VertNormals;
-       errMsg("ntlBlenderDumper","mpTrafo : "<<(*mpTrafo) );
+       //errMsg("ntlBlenderDumper","mpTrafo : "<<(*mpTrafo) );
 
        /* init geometry array, first all standard objects */
        int idCnt = 0;          // give IDs to objects
@@ -107,13 +107,13 @@ int ntlBlenderDumper::renderScene( void )
                if(tid & GEOCLASSTID_SHADER) {
                        ntlGeometryShader *geoshad = (ntlGeometryShader*)(*iter); //dynamic_cast<ntlGeometryShader*>(*iter);
                        hideObjs.push_back( (*iter)->getName() );
-                       geoshad->notifyShaderOfDump(glob->getAniCount(),nrStr,glob->getOutFilename());
+                       geoshad->notifyShaderOfDump(DUMP_FULLGEOMETRY, glob->getAniCount(),nrStr,glob->getOutFilename());
                        for (vector<ntlGeometryObject*>::iterator siter = geoshad->getObjectsBegin();
                                        siter != geoshad->getObjectsEnd();
                                        siter++) {
                                if(debugOut) debMsgStd("ntlBlenderDumper::BuildScene",DM_MSG,"added shader geometry "<<(*siter)->getName(), 8);
 
-                               (*siter)->notifyOfDump(glob->getAniCount(),nrStr,glob->getOutFilename());
+                               (*siter)->notifyOfDump(DUMP_FULLGEOMETRY, glob->getAniCount(),nrStr,glob->getOutFilename(), this->mSimulationTime);
                                bool doDump = false;
                                bool isPreview = false;
                                // only dump final&preview surface meshes
index 7bbae5d..5ab13a9 100644 (file)
@@ -157,9 +157,9 @@ void ntlGeometryObject::initialize(ntlRenderGlobals *glob)
 
 /*! notify object that dump is in progress (e.g. for particles) */
 // default action - do nothing...
-void ntlGeometryObject::notifyOfDump(int frameNr,char *frameNrStr,string outfilename) {
+void ntlGeometryObject::notifyOfDump(int dumtp, int frameNr,char *frameNrStr,string outfilename, double simtime) {
   bool debugOut=false;
-       if(debugOut) debMsgStd("ntlGeometryObject::notifyOfDump",DM_MSG,"obj:"<<this->getName()<<" frame:"<<frameNrStr<<","<<frameNr<<" to "<<outfilename, 10); // DEBUG
+       if(debugOut) debMsgStd("ntlGeometryObject::notifyOfDump",DM_MSG," dt:"<<dumtp<<" obj:"<<this->getName()<<" frame:"<<frameNrStr<<","<<frameNr<<",t"<<simtime<<" to "<<outfilename, 10); // DEBUG
 }
 
 /*****************************************************************************/
index 157d160..ae81d74 100644 (file)
@@ -16,6 +16,8 @@
 class ntlRenderGlobals;
 class ntlTriangle;
 
+#define DUMP_FULLGEOMETRY 1
+#define DUMP_PARTIAL      2
 
 class ntlGeometryObject : public ntlGeometryClass
 {
@@ -38,7 +40,7 @@ class ntlGeometryObject : public ntlGeometryClass
                                vector<ntlVec3Gfx> *normals, int objectId ) = 0;
                
                /*! notify object that dump is in progress (e.g. for particles) */
-               virtual void notifyOfDump(int frameNr,char *frameNrStr,string outfilename);
+               virtual void notifyOfDump(int dumptype, int frameNr,char *frameNrStr,string outfilename, double simtime);
 
                /*! Search the material for this object from the material list */
                void searchMaterial(vector<ntlMaterial *> *mat);
index 7d767b0..6adc562 100644 (file)
@@ -40,7 +40,7 @@ class ntlGeometryShader :
                virtual vector<ntlGeometryObject *>::iterator getObjectsEnd() { return mObjects.end(); }
                
                /*! notify object that dump is in progress (e.g. for field dump) */
-               virtual void notifyShaderOfDump(int frameNr,char *frameNrStr,string outfilename) = 0;
+               virtual void notifyShaderOfDump(int dumptype, int frameNr,char *frameNrStr,string outfilename) = 0;
 
        protected:
 
index 2256e85..131f182 100644 (file)
@@ -380,9 +380,6 @@ private:
        bool mSingleFrameMode;
        //! filename for single frame mode
        string mSingleFrameFilename;
-
-       /*! Two random number streams for photon generation (one for the directions, the other for russion roulette) */
-       //ntlRandomStream *mpRndDirections, *mpRndRoulette;
 };
 
 
index 0a5af90..8286746 100644 (file)
@@ -20,6 +20,8 @@
 #include <zlib.h>
 
 
+// particle object id counter
+int ParticleObjectIdCnt = 1;
 
 /******************************************************************************
  * Standard constructor
@@ -34,19 +36,21 @@ ParticleTracer::ParticleTracer() :
        mPartScale(1.0) , mPartHeadDist( 0.5 ), mPartTailDist( -4.5 ), mPartSegments( 4 ),
        mValueScale(0),
        mValueCutoffTop(0.0), mValueCutoffBottom(0.0),
-       mDumpParts(0), mShowOnly(0), mpTrafo(NULL)
+       mDumpParts(0), mDumpText(0), mDumpTextFile(""), mShowOnly(0), 
+       mNumInitialParts(0), mpTrafo(NULL)
 {
        // check env var
-#ifdef ELBEEM_PLUGIN
-       mDumpParts=1; // default on
-#endif // ELBEEM_PLUGIN
-       if(getenv("ELBEEM_DUMPPARTICLE")) { // DEBUG!
-               int set = atoi( getenv("ELBEEM_DUMPPARTICLE") );
-               if((set>=0)&&(set!=mDumpParts)) {
-                       mDumpParts=set;
-                       debMsgStd("ParticleTracer::notifyOfDump",DM_NOTIFY,"Using envvar ELBEEM_DUMPPARTICLE to set mDumpParts to "<<set<<","<<mDumpParts,8);
-               }
-       }
+//#ifdef ELBEEM_PLUGIN
+       //mDumpParts=1; // default on
+//#endif // ELBEEM_PLUGIN
+       //if(getenv("ELBEEM_DUMPPARTICLE")) { // DEBUG!
+               //int set = atoi( getenv("ELBEEM_DUMPPARTICLE") );
+               //if((set>=0)&&(set!=mDumpParts)) {
+                       //mDumpParts=set;
+                       //debMsgStd("ParticleTrace",DM_NOTIFY,"Using envvar ELBEEM_DUMPPARTICLE to set mDumpParts to "<<set<<","<<mDumpParts,8);
+               //}
+       //}
+       debMsgStd("ParticleTracer::ParticleTracer",DM_MSG,"inited",10);
 };
 
 ParticleTracer::~ParticleTracer() {
@@ -60,13 +64,9 @@ void ParticleTracer::parseAttrList(AttributeList *att)
 {
        AttributeList *tempAtt = mpAttrs; 
        mpAttrs = att;
-       int mNumParticles =0; // UNUSED
-       int mTrailLength  = 0; // UNUSED
-       int mTrailInterval= 0; // UNUSED
-       mNumParticles = mpAttrs->readInt("particles",mNumParticles, "ParticleTracer","mNumParticles", false);
-       mTrailLength  = mpAttrs->readInt("traillength",mTrailLength, "ParticleTracer","mTrailLength", false);
-       mTrailInterval= mpAttrs->readInt("trailinterval",mTrailInterval, "ParticleTracer","mTrailInterval", false);
 
+       mNumInitialParts = mpAttrs->readInt("particles",mNumInitialParts, "ParticleTracer","mNumInitialParts", false);
+       errMsg(" NUMP"," "<<mNumInitialParts);
        mPartScale    = mpAttrs->readFloat("part_scale",mPartScale, "ParticleTracer","mPartScale", false);
        mPartHeadDist = mpAttrs->readFloat("part_headdist",mPartHeadDist, "ParticleTracer","mPartHeadDist", false);
        mPartTailDist = mpAttrs->readFloat("part_taildist",mPartTailDist, "ParticleTracer","mPartTailDist", false);
@@ -75,18 +75,23 @@ void ParticleTracer::parseAttrList(AttributeList *att)
        mValueCutoffTop = mpAttrs->readFloat("part_valcutofftop",mValueCutoffTop, "ParticleTracer","mValueCutoffTop", false);
        mValueCutoffBottom = mpAttrs->readFloat("part_valcutoffbottom",mValueCutoffBottom, "ParticleTracer","mValueCutoffBottom", false);
 
-       mDumpParts    = mpAttrs->readInt  ("part_dump",mDumpParts, "ParticleTracer","mDumpParts", false);
+       mDumpParts   = mpAttrs->readInt  ("part_dump",mDumpParts, "ParticleTracer","mDumpParts", false);
+       mDumpText    = mpAttrs->readInt  ("part_textdump",mDumpText, "ParticleTracer","mDumpText", false);
        mShowOnly    = mpAttrs->readInt  ("part_showonly",mShowOnly, "ParticleTracer","mShowOnly", false);
+       mDumpTextFile= mpAttrs->readString("part_textdumpfile",mDumpTextFile, "ParticleTracer","mDumpTextFile", false);
 
        string matPart;
        matPart = mpAttrs->readString("material_part", "default", "ParticleTracer","material", false);
        setMaterialName( matPart );
-       // trail length has to be at least one, if anything should be displayed
-       //if((mNumParticles>0)&&(mTrailLength<2)) mTrailLength = 2;
+
+       // unused...
+       int mTrailLength  = 0; // UNUSED
+       int mTrailInterval= 0; // UNUSED
+       mTrailLength  = mpAttrs->readInt("traillength",mTrailLength, "ParticleTracer","mTrailLength", false);
+       mTrailInterval= mpAttrs->readInt("trailinterval",mTrailInterval, "ParticleTracer","mTrailInterval", false);
 
        // restore old list
        mpAttrs = tempAtt;
-       //mParts.resize(mTrailLength*mTrailInterval);
 }
 
 /******************************************************************************
@@ -146,7 +151,7 @@ void ParticleTracer::addParticle(float x, float y, float z)
 void ParticleTracer::cleanup() {
        // cleanup
        int last = (int)mParts.size()-1;
-       //for(vector<ParticleObject>::iterator pit= getParticlesBegin();pit!= getParticlesEnd(); pit++) {
+       if(mDumpText>0) { errMsg("ParticleTracer::cleanup","Skipping cleanup due to text dump..."); return; }
 
        for(int i=0; i<=last; i++) {
                if( mParts[i].getActive()==false ) {
@@ -158,47 +163,18 @@ void ParticleTracer::cleanup() {
 }
                
 /******************************************************************************
- * save particle positions before adding a new timestep
- * copy "one index up", newest has to remain unmodified, it will be
- * advanced after the next smiulation step
+ *! dump particles if desired 
  *****************************************************************************/
-void ParticleTracer::savePreviousPositions()
-{
-       //debugOut(" PARTS SIZE "<<mParts.size() ,10);
-       /*
-       if(mTrailIntervalCounter==0) {
-       //errMsg("spp"," PARTS SIZE "<<mParts.size() );
-               for(size_t l=mParts.size()-1; l>0; l--) {
-                       if( mParts[l].size() != mParts[l-1].size() ) {
-                               errFatal("ParticleTracer::savePreviousPositions","Invalid array sizes ["<<l<<"]="<<mParts[l].size()<<
-                                               " ["<<(l+1)<<"]="<<mParts[l+1].size() <<" , total "<< mParts.size() , SIMWORLD_GENERICERROR);
-                               return;
-                       }
-
-                       for(size_t i=0; i<mParts[l].size(); i++) {
-                               mParts[l][i] = mParts[l-1][i];
-                       }
-
-               }
-       } 
-       mTrailIntervalCounter++;
-       if(mTrailIntervalCounter>=mTrailInterval) mTrailIntervalCounter = 0;
-       UNUSED!? */
-}
-
-
-/*! dump particles if desired */
-void ParticleTracer::notifyOfDump(int frameNr,char *frameNrStr,string outfilename) {
+void ParticleTracer::notifyOfDump(int dumptype, int frameNr,char *frameNrStr,string outfilename, double simtime) {
        debMsgStd("ParticleTracer::notifyOfDump",DM_MSG,"obj:"<<this->getName()<<" frame:"<<frameNrStr, 10); // DEBUG
 
-       if(mDumpParts>0) {
+       if(
+                       (dumptype==DUMP_FULLGEOMETRY)&&
+                       (mDumpParts>0)) {
                // dump to binary file
                std::ostringstream boutfilename("");
-               //boutfilename << ecrpath.str() << glob->getOutFilename() <<"_"<< this->getName() <<"_" << frameNrStr << ".obj";
-               //boutfilename << outfilename <<"_particles_"<< this->getName() <<"_" << frameNrStr<< ".gz";
                boutfilename << outfilename <<"_particles_" << frameNrStr<< ".gz";
-               debMsgStd("ParticleTracer::notifyOfDump",DM_MSG,"B-Dumping: "<< this->getName() 
-                               <<", particles:"<<mParts.size()<<" "<< " to "<<boutfilename.str()<<" #"<<frameNr , 7);
+               debMsgStd("ParticleTracer::notifyOfDump",DM_MSG,"B-Dumping: "<< this->getName() <<", particles:"<<mParts.size()<<" "<< " to "<<boutfilename.str()<<" #"<<frameNr , 7);
 
                // output to zipped file
                gzFile gzf;
@@ -210,21 +186,26 @@ void ParticleTracer::notifyOfDump(int frameNr,char *frameNrStr,string outfilenam
                        numParts = 0;
                        for(size_t i=0; i<mParts.size(); i++) {
                                if(!mParts[i].getActive()) continue;
+                               //if(mParts[i].getLifeTime()<30) { continue; } //? CHECK
                                numParts++;
                        }
                        gzwrite(gzf, &numParts, sizeof(numParts));
                        for(size_t i=0; i<mParts.size(); i++) {
                                if(!mParts[i].getActive()) { continue; }
-                               if(mParts[i].getLifeTime()<30) { continue; } //? CHECK
+                               //if(mParts[i].getLifeTime()<30) { continue; } //? CHECK
                                ParticleObject *p = &mParts[i];
-                               int type = p->getType();
+                               //int type = p->getType();  // export whole type info
+                               int type = p->getFlags(); // debug export whole type & status info
                                ntlVec3Gfx pos = p->getPos();
                                float size = p->getSize();
 
                                if(type&PART_FLOAT) { // WARNING same handling for dump!
                                        // add one gridcell offset
                                        //pos[2] += 1.0; 
-                               }
+                               } 
+                               // display as drop for now externally
+                               //else if(type&PART_TRACER) { type |= PART_DROP; }
+
                                pos = (*mpTrafo) * pos;
 
                                ntlVec3Gfx v = p->getVel();
@@ -240,6 +221,64 @@ void ParticleTracer::notifyOfDump(int frameNr,char *frameNrStr,string outfilenam
                        gzclose( gzf );
                }
        } // dump?
+
+       // dfor partial & full dump
+       if(mDumpText>0) {
+               // dump to binary file
+               std::ostringstream boutfilename("");
+               if(mDumpTextFile.length()>1) {   
+                       boutfilename << mDumpTextFile <<  ".cpart2"; 
+               } else {                           
+                       boutfilename << outfilename <<"_particles" <<  ".cpart2"; 
+               }
+               debMsgStd("ParticleTracer::notifyOfDump",DM_MSG,"T-Dumping: "<< this->getName() <<", particles:"<<mParts.size()<<" "<< " to "<<boutfilename.str()<<" #"<<frameNr , 7);
+
+               int numParts = 0;
+               // only dump bubble particles
+               for(size_t i=0; i<mParts.size(); i++) {
+                       //if(!mParts[i].getActive()) continue;
+                       //if(!(mParts[i].getType()&PART_BUBBLE)) continue;
+                       numParts++;
+               }
+
+               // output to text file
+               gzFile gzf;
+               if(frameNr==0) {
+                       gzf = gzopen(boutfilename.str().c_str(), "w0");
+
+                       gzprintf( gzf, "\n\n# cparts generated by elbeem \n# no. of parts \nN %d \n\n",numParts);
+                       // fixed time scale for now
+                       gzprintf( gzf, "T %f \n\n", 1.0);
+               } else {
+                       gzf = gzopen(boutfilename.str().c_str(), "a+0");
+               }
+
+               // add current set
+               if(gzf) {
+                       gzprintf( gzf, "\n\n# new set at frame %d,t%f,p%d --------------------------------- \n\n", frameNr, simtime, numParts );
+                       gzprintf( gzf, "S %f \n\n", simtime );
+                       
+                       for(size_t i=0; i<mParts.size(); i++) {
+                               ParticleObject *p = &mParts[i];
+                               ntlVec3Gfx pos = p->getPos();
+                               float size = p->getSize();
+                               if(!mParts[i].getActive()) { size=0.; } // switch "off"
+
+                               pos = (*mpTrafo) * pos;
+                               ntlVec3Gfx v = p->getVel();
+                               v[0] *= mpTrafo->value[0][0];
+                               v[1] *= mpTrafo->value[1][1];
+                               v[2] *= mpTrafo->value[2][2];
+                               
+                               gzprintf( gzf, "P %f %f %f \n", pos[0],pos[1],pos[2] );
+                               gzprintf( gzf, "s %f \n", size );
+                               gzprintf( gzf, "\n", size );
+                       }
+
+                       gzprintf( gzf, "# %d end  ", frameNr );
+                       gzclose( gzf );
+               }
+       }
 }
 
 
@@ -293,6 +332,7 @@ void ParticleTracer::getTriangles( vector<ntlTriangle> *triangles,
                        case 2: if(!(type&PART_DROP))   continue; break;
                        case 3: if(!(type&PART_INTER))  continue; break;
                        case 4: if(!(type&PART_FLOAT))  continue; break;
+                       case 5: if(!(type&PART_TRACER))  continue; break;
                        }
                } else {
                        // by default dont display inter
index b9179b6..453c1a7 100644 (file)
@@ -16,24 +16,31 @@ template<class Scalar> class ntlMatrix4x4;
 #define PART_DROP   (1<< 2)
 #define PART_INTER  (1<< 3)
 #define PART_FLOAT  (1<< 4)
+#define PART_TRACER (1<< 5)
 
 // particle state
 #define PART_IN     (1<< 8)
 #define PART_OUT    (1<< 9)
 #define PART_INACTIVE (1<<10)
 
+// defines for particle movement
+#define MOVE_FLOATS 1
+#define FLOAT_JITTER 0.03
+
+extern int ParticleObjectIdCnt;
+
 //! A single particle
 class ParticleObject
 {
        public:
        //! Standard constructor
        inline ParticleObject(ntlVec3Gfx mp) :
-                       mPos(mp),mVel(0.0), mSize(1.0), mStatus(0),mLifeTime(0) { };
+                       mPos(mp),mVel(0.0), mSize(1.0), mStatus(0),mLifeTime(0) { mId = ParticleObjectIdCnt++; };
        //! Copy constructor
        inline ParticleObject(const ParticleObject &a) :
                        mPos(a.mPos), mVel(a.mVel), mSize(a.mSize), 
                        mStatus(a.mStatus),
-                       mLifeTime(a.mLifeTime) { };
+                       mLifeTime(a.mLifeTime) { mId = ParticleObjectIdCnt++; };
        //! Destructor
        inline ~ParticleObject() { /* empty */ };
 
@@ -81,8 +88,12 @@ class ParticleObject
                //! set type (lower byte)
                inline void setLifeTime(int set) { mLifeTime = set; }
                
+               inline int getId() const { return mId; }
+
        protected:
 
+               /*! only for debugging */
+               int mId;
                /*! the particle position */
                ntlVec3Gfx mPos;
                /*! the particle velocity */
@@ -109,9 +120,6 @@ class ParticleTracer :
                //! add a particle at this position
                void addParticle(float x, float y, float z);
 
-               //! save particle positions before adding a new timestep
-               void savePreviousPositions();
-
                //! draw the particle array
                void draw();
                
@@ -125,6 +133,9 @@ class ParticleTracer :
                
                //! get the number of particles
                inline int  getNumParticles()                           { return mParts.size(); }
+               //! set/get the number of particles
+               inline void setNumInitialParticles(int set) { mNumInitialParts=set; }
+               inline int  getNumInitialParticles()          { return mNumInitialParts; }
 
                //! iterate over all newest particles (for advancing positions)
                inline vector<ParticleObject>::iterator getParticlesBegin() { return mParts.begin(); }
@@ -149,7 +160,13 @@ class ParticleTracer :
                
                /*! set/get dump flag */
                inline void setDumpParts(bool set) { mDumpParts = set; }
-               inline bool getDumpParts() { return mDumpParts; }
+               inline bool getDumpParts()         { return mDumpParts; }
+               /*! set/get dump flag */
+               inline void setDumpText(bool set) { mDumpText = set; }
+               inline bool getDumpText()         { return mDumpText; }
+               /*! set/get dump text file */
+               inline void setDumpTextFile(std::string set) { mDumpTextFile = set; }
+               inline std::string getDumpTextFile()         { return mDumpTextFile; }
                
                //! set the particle scaling factor
                inline void setPartScale(float set) { mPartScale = set; }
@@ -161,7 +178,7 @@ class ParticleTracer :
                                vector<ntlVec3Gfx> *vertices, 
                                vector<ntlVec3Gfx> *normals, int objectId );
 
-               virtual void notifyOfDump(int frameNr,char *frameNrStr,string outfilename);
+               virtual void notifyOfDump(int dumptype, int frameNr,char *frameNrStr,string outfilename,double simtime);
                // free deleted particles
                void cleanup();
 
@@ -194,8 +211,14 @@ class ParticleTracer :
 
                /*! dump particles (or certain types of) to disk? */
                int mDumpParts;
+               /*! dump particles (or certain types of) to disk? */
+               int mDumpText;
+               /*! text dump output file */
+               std::string mDumpTextFile;
                /*! show only a certain type (debugging) */
                int mShowOnly;
+               /*! no. of particles to init */
+               int mNumInitialParts;
 
                //! transform matrix
                ntlMatrix4x4<gfxReal> *mpTrafo;
index 481eb83..116d68d 100644 (file)
@@ -156,6 +156,7 @@ int SimulationObject::initializeLbmSimulation(ntlRenderGlobals *glob)
        mpLbm->setGeoEnd( mGeoEnd );
        mpLbm->setRenderGlobals( mpGlob );
        mpLbm->setName( getName() + "_lbm" );
+       mpLbm->setParticleTracer( mpParts );
        if(mpElbeemSettings) {
                // set further settings from API struct init
                mpLbm->setSmoothing(1.0 * mpElbeemSettings->surfaceSmoothing, 1.0 * mpElbeemSettings->surfaceSmoothing);
@@ -173,6 +174,7 @@ int SimulationObject::initializeLbmSimulation(ntlRenderGlobals *glob)
                mpLbm->setDomainBound(dinitType);
                mpLbm->setDomainPartSlip(mpElbeemSettings->obstaclePartslip);
                mpLbm->setDumpVelocities(mpElbeemSettings->generateVertexVectors);
+               mpLbm->setFarFieldSize(mpElbeemSettings->farFieldSize);
                debMsgStd("SimulationObject::initialize",DM_MSG,"Added domain bound: "<<dinitType<<" ps="<<mpElbeemSettings->obstaclePartslip<<" vv"<<mpElbeemSettings->generateVertexVectors<<","<<mpLbm->getDumpVelocities(), 9 );
 
                debMsgStd("SimulationObject::initialize",DM_MSG,"Set ElbeemSettings values "<<mpLbm->getGenerateParticles(),10);
@@ -236,7 +238,6 @@ int SimulationObject::initializeLbmSimulation(ntlRenderGlobals *glob)
        mpParts->setCastShadows( false );
        mpParts->setReceiveShadows( false );
        mpParts->searchMaterial( glob->getMaterials() );
-       mpLbm->initParticles(mpParts);
 
        // this has to be inited here - before, the values might be unknown
        ntlGeometryObject *surf = mpLbm->getSurfaceGeoObj();
@@ -288,9 +289,6 @@ void SimulationObject::step( void )
        if(mpParam->getCurrentAniFrameTime()>0.0) {
                // dont advance for stopped time
                mpLbm->step();
-
-               mpParts->savePreviousPositions();
-               mpLbm->advanceParticles(mpParts);
                mTime += mpParam->getTimestep();
        }
        if(mpLbm->getPanic()) mPanic = true;
@@ -399,8 +397,8 @@ void SimulationObject::setMouseClick()
 }
 
 /*! notify object that dump is in progress (e.g. for field dump) */
-void SimulationObject::notifyShaderOfDump(int frameNr,char *frameNrStr,string outfilename) {
+void SimulationObject::notifyShaderOfDump(int dumptype, int frameNr,char *frameNrStr,string outfilename) {
        if(!mpLbm) return;
-       mpLbm->notifySolverOfDump(frameNr,frameNrStr,outfilename);
+       mpLbm->notifySolverOfDump(dumptype, frameNr,frameNrStr,outfilename);
 }
 
index 4e7bc19..b8d1e90 100644 (file)
@@ -91,7 +91,7 @@ class SimulationObject :
                virtual int postGeoConstrInit(ntlRenderGlobals *glob) { return initializeLbmSimulation(glob); };
                virtual int initializeShader() { /* ... */ return true; };
                /*! notify object that dump is in progress (e.g. for field dump) */
-               virtual void notifyShaderOfDump(int frameNr,char *frameNrStr,string outfilename);
+               virtual void notifyShaderOfDump(int dumptype, int frameNr,char *frameNrStr,string outfilename);
                /*! simluation interface: draw the simulation with OpenGL */
                virtual void draw( void ) {};
                virtual vector<ntlGeometryObject *>::iterator getObjectsBegin();
index 7217964..b241834 100644 (file)
@@ -211,7 +211,7 @@ class LbmFsgrSolver :
                //! finish the init with config file values (allocate arrays...) 
                virtual bool initializeSolver(); 
                //! notify object that dump is in progress (e.g. for field dump) 
-               virtual void notifySolverOfDump(int frameNr,char *frameNrStr,string outfilename);
+               virtual void notifySolverOfDump(int dumptype, int frameNr,char *frameNrStr,string outfilename);
 
 #if LBM_USE_GUI==1
                //! show simulation info (implement LbmSolverInterface pure virtual func)
@@ -267,9 +267,9 @@ class LbmFsgrSolver :
                /* simulation object interface, just calls stepMain */
                virtual void step();
                /*! init particle positions */
-               virtual int initParticles(ParticleTracer *partt);
+               virtual int initParticles();
                /*! move all particles */
-               virtual void advanceParticles(ParticleTracer *partt );
+               virtual void advanceParticles();
 
 
                /*! debug object display (used e.g. for preview surface) */
@@ -436,8 +436,6 @@ class LbmFsgrSolver :
                int mForceTadapRefine;
                //! border cutoff value
                int mCutoff;
-               //! store particle tracer
-               ParticleTracer *mpParticles;
 
                // strict debug interface
 #              if FSGR_STRICT_DEBUG==1
@@ -460,11 +458,10 @@ class LbmFsgrSolver :
                void initTestdata();
                void destroyTestdata();
                void handleTestdata();
-               void exportTestdata();
                void set3dHeight(int ,int );
        public:
                // needed from testdata
-               void find3dHeight(int i,int j, LbmFloat prev, LbmFloat &ret, LbmFloat &retux, LbmFloat &retuy);
+               void find3dHeight(int i,int j, LbmFloat prev, LbmFloat &ret, LbmFloat *retux, LbmFloat *retuy);
 #endif // LBM_INCLUDE_TESTSOLVERS==1
 
        public: // former LbmModelLBGK  functions
@@ -700,9 +697,8 @@ class LbmFsgrSolver :
 #define _RFLAG_NB(level,xx,yy,zz,set, dir) mLevel[level].mprsFlags[set][ LBMGI((level),(xx)+this->dfVecX[dir],(yy)+this->dfVecY[dir],(zz)+this->dfVecZ[dir],set) ]
 #define _RFLAG_NBINV(level,xx,yy,zz,set, dir) mLevel[level].mprsFlags[set][ LBMGI((level),(xx)+this->dfVecX[this->dfInv[dir]],(yy)+this->dfVecY[this->dfInv[dir]],(zz)+this->dfVecZ[this->dfInv[dir]],set) ]
 
-// array data layouts
-// standard array layout  -----------------------------------------------------------------------------------------------
-#define ALSTRING "Standard Array Layout"
+// array handling  -----------------------------------------------------------------------------------------------
+
 //#define _LBMQI(level, ii,ij,ik, is, lunused) ( ((is)*mLevel[level].lOffsz) + (mLevel[level].lOffsy*(ik)) + (mLevel[level].lOffsx*(ij)) + (ii) )
 #define _LBMQI(level, ii,ij,ik, is, lunused) ( (mLevel[level].lOffsy*(ik)) + (mLevel[level].lOffsx*(ij)) + (ii) )
 #define _QCELL(level,xx,yy,zz,set,l) (mLevel[level].mprsCells[(set)][ LBMQI((level),(xx),(yy),(zz),(set), l)*dTotalNum +(l)])
index abe9ef6..a417c33 100644 (file)
@@ -396,8 +396,6 @@ LbmFsgrSolver::LbmFsgrSolver() :
                mGaussw[n] = mGaussw[n]/totGaussw;
        }
 
-       mpParticles = NULL;
-       //addDrop(false,0,0);
 }
 
 /*****************************************************************************/
@@ -486,6 +484,8 @@ void LbmFsgrSolver::parseAttrList()
 #else // LBM_INCLUDE_TESTSOLVERS!=1
        // off by default
        mUseTestdata = 0;
+       if(mFarFieldSize>=2.) mUseTestdata=1; // equiv. to test solver check
+       if(mUseTestdata) { mMaxRefine=0; } // force fsgr off
 #endif // LBM_INCLUDE_TESTSOLVERS!=1
 }
 
@@ -573,7 +573,14 @@ void LbmFsgrSolver::initLevelOmegas()
  *****************************************************************************/
 bool LbmFsgrSolver::initializeSolver()
 {
-//  debMsgStd("LbmFsgrSolver::initialize",DM_MSG,"Init start... (Layout:"<<ALSTRING<<") "<<this->mInitDone<<" "<<((int)this),1);
+  debMsgStd("LbmFsgrSolver::initialize",DM_MSG,"Init start... "<<this->mInitDone<<" "<<(void*)this,1);
+
+       // init cppf stage
+       if(mCppfStage>0) {
+               this->mSizex *= mCppfStage;
+               this->mSizey *= mCppfStage;
+               this->mSizez *= mCppfStage;
+       }
 
        // size inits to force cubic cells and mult4 level dimensions
        // and make sure we dont allocate too much...
@@ -584,9 +591,6 @@ bool LbmFsgrSolver::initializeSolver()
        double sizeReduction = 1.0;
        double memEstFromFunc = -1.0;
        string memreqStr("");   
-#if LBM_INCLUDE_TESTSOLVERS==1
-       if(mUseTestdata) { mMaxRefine=0; } // force fsgr off
-#endif
        while(!memOk) {
                initGridSizes( this->mSizex, this->mSizey, this->mSizez,
                                this->mvGeoStart, this->mvGeoEnd, mMaxRefine, PARALLEL);
@@ -750,16 +754,25 @@ bool LbmFsgrSolver::initializeSolver()
        mMaxTimestep = this->mpParam->getTimestep();
 
        // init isosurf
+       this->mpIso->setIsolevel( this->mIsoValue );
 #if LBM_INCLUDE_TESTSOLVERS==1
        if(mUseTestdata) {
                mpTest->setMaterialName( this->mpIso->getMaterialName() );
                delete this->mpIso;
                this->mpIso = mpTest;
+               if(mpTest->mDebugvalue1>0.0) { // 3d off
+                       mpTest->setIsolevel(-100.0);
+               } else {
+                       mpTest->setIsolevel( this->mIsoValue );
+               }
        }
 #endif // ELBEEM_PLUGIN!=1
-       this->mpIso->setIsolevel( this->mIsoValue );
        // approximate feature size with mesh resolution
        float featureSize = mLevel[ mMaxRefine ].nodeSize*0.5;
+       // smooth vars defined in solver_interface, set by simulation object
+       // reset for invalid values...
+       if((this->mSmoothSurface<0.)||(this->mSmoothSurface>50.)) this->mSmoothSurface = 1.;
+       if((this->mSmoothNormals<0.)||(this->mSmoothNormals>50.)) this->mSmoothNormals = 1.;
        this->mpIso->setSmoothSurface( this->mSmoothSurface * featureSize );
        this->mpIso->setSmoothNormals( this->mSmoothNormals * featureSize );
 
@@ -984,6 +997,12 @@ bool LbmFsgrSolver::initializeSolver()
                        errMsg("LbmFsgrSolver::init","No preview in 2D allowed!");
                        this->mOutputSurfacePreview = 0; }
        }
+#if LBM_USE_GUI==1
+       if(this->mOutputSurfacePreview) {
+               errMsg("LbmFsgrSolver::init","No preview in GUI mode... mOutputSurfacePreview=0");
+               this->mOutputSurfacePreview = 0; }
+#endif // LBM_USE_GUI==1
+       
        if(this->mOutputSurfacePreview) {
 
                // same as normal one, but use reduced size
@@ -1017,13 +1036,17 @@ bool LbmFsgrSolver::initializeSolver()
 
 
        // now really done...
-  debugOut("LbmFsgrSolver::initialize : Init done ...",10);
+  debugOut("LbmFsgrSolver::initialize : SurfaceGen: SmsOrg("<<this->mSmoothSurface<<","<<this->mSmoothNormals<<","<<featureSize<<"), Iso("<<this->mpIso->getSmoothSurface()<<","<<this->mpIso->getSmoothNormals()<<") ",10);
+  debugOut("LbmFsgrSolver::initialize : Init done ... ",10);
        this->mInitDone = 1;
 
 #if LBM_INCLUDE_TESTSOLVERS==1
        initTestdata();
 #endif // ELBEEM_PLUGIN!=1
+       // not inited? dont use...
+       if(mCutoff<0) mCutoff=0;
 
+       initParticles();
        return true;
 }
 
index 04e0456..ee2f8f2 100644 (file)
@@ -34,7 +34,7 @@ LbmSolverInterface::LbmSolverInterface() :
        mBoundarySouth( (CellFlagType)(CFBnd) ),mBoundaryTop(  (CellFlagType)(CFBnd) ),mBoundaryBottom( (CellFlagType)(CFBnd) ),
   mInitDone( false ),
   mInitDensityGradient( false ),
-       mpAttrs( NULL ), mpParam( NULL ),
+       mpAttrs( NULL ), mpParam( NULL ), mpParticles(NULL),
        mNumParticlesLost(0), 
        mNumInvalidDfs(0), mNumFilledCells(0), mNumEmptiedCells(0), mNumUsedCells(0), mMLSUPS(0),
        mDebugVelScale( 0.01 ), mNodeInfoString("+"),
@@ -53,7 +53,7 @@ LbmSolverInterface::LbmSolverInterface() :
        mDumpVelocities(false),
        mMarkedCells(), mMarkedCellIndex(0),
        mDomainBound("noslip"), mDomainPartSlipValue(0.1),
-       mTForceStrength(0.0)
+       mTForceStrength(0.0), mFarFieldSize(0.), mCppfStage(0)
 {
 #if ELBEEM_PLUGIN==1
        if(gDebugLevel<=1) mSilent = true;
@@ -242,7 +242,12 @@ void LbmSolverInterface::parseStdAttrList() {
        
        // new test vars
        mTForceStrength = mpAttrs->readFloat("tforcestrength", mTForceStrength,"LbmSolverInterface", "mTForceStrength", false);
+       mFarFieldSize = mpAttrs->readFloat("farfieldsize", mFarFieldSize,"LbmSolverInterface", "mFarFieldSize", false);
+       // old compat
+       float sizeScale = mpAttrs->readFloat("test_scale", 0.,"LbmTestdata", "mSizeScale", false);
+       if((mFarFieldSize<=0.)&&(sizeScale>0.)) { mFarFieldSize=sizeScale; errMsg("LbmTestdata","Warning - using mSizeScale..."); }
 
+       mCppfStage = mpAttrs->readFloat("cppfstage", mCppfStage,"LbmSolverInterface", "mCppfStage", false);
        mPartGenProb = mpAttrs->readFloat("partgenprob", mPartGenProb,"LbmFsgrSolver", "mPartGenProb", false);
 }
 
index c86bbce..ee4ecd7 100644 (file)
@@ -260,7 +260,7 @@ class LbmSolverInterface
                virtual bool initializeSolver() =0;
                
                /*! notify object that dump is in progress (e.g. for field dump) */
-               virtual void notifySolverOfDump(int frameNr,char *frameNrStr,string outfilename) = 0;
+               virtual void notifySolverOfDump(int dumptype, int frameNr,char *frameNrStr,string outfilename) = 0;
 
                /*! parse a boundary flag string */
                CellFlagType readBoundaryFlagInt(string name, int defaultValue, string source,string target, bool needed);
@@ -273,8 +273,8 @@ class LbmSolverInterface
                virtual void prepareVisualization() { /* by default off */ };
 
                /*! particle handling */
-               virtual int initParticles(ParticleTracer *partt) = 0;
-               virtual void advanceParticles(ParticleTracer *partt ) = 0;
+               virtual int initParticles() = 0;
+               virtual void advanceParticles() = 0;
                /*! get surface object (NULL if no surface) */
                ntlGeometryObject* getSurfaceGeoObj() { return mpIso; }
 
@@ -329,6 +329,9 @@ class LbmSolverInterface
                inline void setParametrizer(Parametrizer *set) { mpParam = set; }
                /*! get parametrizer pointer */
                inline Parametrizer *getParametrizer() { return mpParam; }
+               /*! get/set particle pointer */
+               inline void setParticleTracer(ParticleTracer *set) { mpParticles = set; }
+               inline ParticleTracer *getParticleTracer() { return mpParticles; }
 
                /*! set density gradient init from e.g. init test cases */
                inline void setInitDensityGradient(bool set) { mInitDensityGradient = set; }
@@ -379,6 +382,12 @@ class LbmSolverInterface
                //! set/get dump velocities flag
                inline void setDomainPartSlip(LbmFloat set)     { mDomainPartSlipValue = set; }
                inline LbmFloat getDomainPartSlip() const       { return mDomainPartSlipValue; }
+               //! set/get far field size
+               inline void setFarFieldSize(LbmFloat set)       { mFarFieldSize = set; }
+               inline LbmFloat getFarFieldSize() const { return mFarFieldSize; }
+               //! set/get cp stage
+               inline void setCpStage(int set) { mCppfStage = set; }
+               inline int getCpStage() const   { return mCppfStage; }
 
                // cell iterator interface
                
@@ -474,6 +483,8 @@ class LbmSolverInterface
 
                /*! get parameters from this parametrize in finishInit */
                Parametrizer *mpParam;
+               //! store particle tracer
+               ParticleTracer *mpParticles;
 
                /*! number of particles lost so far */
                int mNumParticlesLost;
@@ -553,6 +564,10 @@ class LbmSolverInterface
                //! test vars
                // strength of applied force
                LbmFloat mTForceStrength;
+               // 
+               LbmFloat mFarFieldSize;
+               // 
+               int mCppfStage;
 
 };
 
index c8f0678..10fca60 100644 (file)
@@ -108,7 +108,7 @@ void LbmFsgrSolver::stepMain()
                int avgcls = (int)(mAvgNumUsedCells/(LONGINT)this->mStepCnt);
        debMsgStd("LbmFsgrSolver::step", DM_MSG, this->mName<<" cnt:"<<this->mStepCnt<<" t:"<<mSimulationTime<<
                //debMsgDirect( 
-                       "mlsups(curr:"<<this->mMLSUPS<<
+                       " mlsups(curr:"<<this->mMLSUPS<<
                        " avg:"<<(mAvgMLSUPS/mAvgMLSUPSCnt)<<"), "<< sepStr<<
                        " totcls:"<<(this->mNumUsedCells)<< sepStr<<
                        " avgcls:"<< avgcls<< sepStr<<
@@ -122,9 +122,6 @@ void LbmFsgrSolver::stepMain()
                        " probs:"<<mNumProblems<< sepStr<<
                        " simt:"<<mSimulationTime<< sepStr<<
                        " for '"<<this->mName<<"' " , 10);
-
-               debMsgDirect(std::endl);
-               debMsgDirect(this->mStepCnt<<": dccd="<< mCurrentMass<<"/"<<mCurrentVolume<<"(fix="<<this->mFixMass<<",ini="<<mInitialMass<<") ");
                debMsgDirect(std::endl);
 
                // nicer output
@@ -208,6 +205,9 @@ void LbmFsgrSolver::stepMain()
        if((mUseTestdata)&&(this->mInitDone)) { handleTestdata(); }
 #endif
 
+       // advance positions with current grid
+       advanceParticles();
+
        // one of the last things to do - adapt timestep
        // was in fineAdvance before... 
        if(mTimeAdap) {
@@ -220,6 +220,13 @@ void LbmFsgrSolver::stepMain()
        if( (!finite(mMxvx)) || (!finite(mMxvy)) || (!finite(mMxvz)) ) { CAUSE_PANIC; }
        if( (!finite(mCurrentMass)) || (!finite(mCurrentVolume)) ) { CAUSE_PANIC; }
 #endif // WIN32
+
+       // output total step time
+       timeend = getTime();
+       debMsgStd("LbmFsgrSolver::stepMain",DM_MSG,"step:"<<this->mStepCnt
+                       <<": dccd="<< mCurrentMass<<"/"<<mCurrentVolume<<"(fix="<<this->mFixMass<<",ini="<<mInitialMass<<"), "
+                       <<" totst:"<<getTimeString(timeend-timestart), 3);
+
  //#endif // ELBEEM_PLUGIN!=1
 }
 
@@ -283,11 +290,12 @@ LbmFsgrSolver::mainLoop(int lev)
        int      calcCellsFilled = this->mNumFilledCells;
        int      calcCellsEmptied = this->mNumEmptiedCells;
        int      calcNumUsedCells = this->mNumUsedCells;
+       const int cutMin  = 1;
+       const int cutConst = mCutoff+1;
 
 #      if LBM_INCLUDE_TESTSOLVERS==1
-       if((mUseTestdata)&&(mpTest->mDebugvalue1>0.0)) { 
-               // 3d region off... quit
-               this->mpIso->setIsolevel(-100.0); return; }
+       // 3d region off... quit
+       if((mUseTestdata)&&(mpTest->mDebugvalue1>0.0)) { return; }
 #endif // ELBEEM_PLUGIN!=1
        //printLbmCell(lev, 6,6,16, mLevel[lev].setCurr ); // DEBUG
        
@@ -814,13 +822,18 @@ LbmFsgrSolver::mainLoop(int lev)
                // for inflow no pargen test
                // NOBUBBB!
                if((this->mInitDone) //&&(mUseTestdata) 
-                               && (!((oldFlag|newFlag)&CFNoNbEmpty)) 
+                               //&& (!((oldFlag|newFlag)&CFNoNbEmpty)) 
                                && (!((oldFlag|newFlag)&CFNoDelete)) 
                                && (this->mPartGenProb>0.0)) {
+                       bool doAdd = true;
                        bool bndOk=true;
-                       if( (i<1+mCutoff)||(i>this->mSizex-1-mCutoff)||
-                                       (j<1+mCutoff)||(j>this->mSizey-1-mCutoff)||
-                                       (k<1+mCutoff)||(k>this->mSizez-1-mCutoff) ) { bndOk=false; }
+                       //if( (i<1+mCutoff)||(i>this->mSizex-1-mCutoff)||
+                                       //(j<1+mCutoff)||(j>this->mSizey-1-mCutoff)||
+                                       //(k<1+mCutoff)||(k>this->mSizez-1-mCutoff) ) { bndOk=false; }
+                       if( (i<cutMin)||(i>this->mSizex-cutMin)||
+                                       (j<cutMin)||(j>this->mSizey-cutMin)||
+                                       (k<cutMin)||(k>this->mSizez-cutMin) ) { bndOk=false; }
+                       if(!bndOk) doAdd=false;
                        
                        LbmFloat realWorldFac = (mLevel[lev].simCellSize / mLevel[lev].timestep);
                        LbmFloat rux = (ux * realWorldFac);
@@ -828,12 +841,19 @@ LbmFsgrSolver::mainLoop(int lev)
                        LbmFloat ruz = (uz * realWorldFac);
                        LbmFloat rl = norm(ntlVec3Gfx(rux,ruy,ruz));
                        // WHMOD
-                       const LbmFloat val2fac = 1.0; //? TODO N test? /(this->mPartGenProb); // FIXME remove factor!
-                       bool doAdd = true;
 
                        LbmFloat prob = (rand()/(RAND_MAX+1.0));
                        LbmFloat basethresh = this->mPartGenProb*lcsmqo*rl;
-                       if( (prob< (basethresh*rl)) && (lcsmqo>0.0095) && (rl>2.5) ) {
+
+                       // reduce probability in outer region
+                       if( (i<cutConst)||(i>this->mSizex-cutConst) ){ prob *= 0.5; }
+                       if( (j<cutConst)||(j>this->mSizey-cutConst) ){ prob *= 0.5; }
+                       if( (k<cutConst)||(k>this->mSizez-cutConst) ){ prob *= 0.5; }
+
+//#define RWVEL_THRESH 1.0
+#define RWVEL_THRESH 1.5
+#define RWVEL_WINDTHRESH (RWVEL_THRESH*0.5)
+                       if( (prob< (basethresh*rl)) && (lcsmqo>0.0095) && (rl>RWVEL_THRESH) ) {
                                // add
                        } else {
                                doAdd = false; // dont...
@@ -851,12 +871,12 @@ LbmFsgrSolver::mainLoop(int lev)
 
                        // "wind" disturbance
                        // use realworld relative velocity here instead?
-                       if( 
-                                       ((rl>1.0) && (lcsmqo<P_LCSMQO)) // normal checks
-                                       ||(k>this->mSizez-SLOWDOWNREGION) ) {
+                       if( (doAdd) && (
+                                       ((rl>RWVEL_WINDTHRESH) && (lcsmqo<P_LCSMQO)) // normal checks
+                                       ||(k>this->mSizez-SLOWDOWNREGION)   )) {
                                LbmFloat nuz = uz;
                                LbmFloat jdf; // = 0.05 * (rand()/(RAND_MAX+1.0));
-                               if(rl>1.0) jdf *= (rl-1.0);
+                               if(rl>RWVEL_WINDTHRESH) jdf *= (rl-RWVEL_WINDTHRESH);
                                if(k>this->mSizez-SLOWDOWNREGION) {
                                        // special case
                                        LbmFloat zfac = (LbmFloat)( k-(this->mSizez-SLOWDOWNREGION) );
@@ -881,7 +901,6 @@ LbmFsgrSolver::mainLoop(int lev)
                        if(usqr<0.0001) doAdd=false;   // TODO check!?
                        // if outside, and 20% above sea level, delete, TODO really check level?
                        //if((!bndOk)&&((LbmFloat)k>pTest->mFluidHeight*1.5)) { doAdd=true; } // FORCEDISSOLVE
-                       //? if(!bndOk) doAdd=false;
                        //if(this->mStepCnt>700) errMsg("DFJITT"," at "<<PRINT_IJK<<"rwl:"<<rl<<"  usqr:"<<usqr <<" qo:"<<lcsmqo<<" add="<<doAdd );
 
                        if( (doAdd)  ) { // ADD DROP
@@ -911,15 +930,14 @@ LbmFsgrSolver::mainLoop(int lev)
                                LbmFloat srcj = j+0.5+jpy; // TEST? + (pv[1]*1.41);
                                LbmFloat srck = k+0.5+jpz; // TEST? + (pv[2]*1.41);
                                int type=0;
-                               //if((s%3)!=2) {
-
-                                       type=PART_DROP;
-                                       // drop
-                                       srci += (pv[0]*1.41);
-                                       srcj += (pv[1]*1.41);
-                                       srck += (pv[2]*1.41);
-                                       if(!(RFLAG(lev, (int)(srci),(int)(srcj),(int)(srck),SRCS(lev)) &CFEmpty)) continue; // only add in good direction
-                               //} else { type=PART_FLOAT; }
+                               //if((s%3)!=2) {} else { type=PART_FLOAT; }
+                               //type = PART_DROP;
+                               type = PART_INTER;
+                               // drop
+                               /*srci += (pv[0]*1.41);
+                               srcj += (pv[1]*1.41);
+                               srck += (pv[2]*1.41);
+                               if(!(RFLAG(lev, (int)(srci),(int)(srcj),(int)(srck),SRCS(lev)) &CFEmpty)) continue; // only add in good direction */
 
                                pv *= len;
                                LbmFloat size = 1.0+ 9.0* (rand()/(RAND_MAX+1.0));
@@ -929,18 +947,20 @@ LbmFsgrSolver::mainLoop(int lev)
                                //? mpParticles->getLast()->advanceVel(); // advance a bit outwards
                                mpParticles->getLast()->setStatus(PART_IN);
                                mpParticles->getLast()->setType(type);
-                               //mpParticles->getLast()->setType(PART_INTER);
                                //if((s%3)==2) mpParticles->getLast()->setType(PART_FLOAT);
                                mpParticles->getLast()->setSize(size);
                                //errMsg("NEWPART"," at "<<PRINT_IJK<<"   u="<<PRINT_VEC(ux,uy,uz) <<" RWu="<<PRINT_VEC(rux,ruy,ruz)<<" add"<<doAdd<<" pvel="<<pv );
                                //errMsg("NEWPART"," at "<<PRINT_IJK<<"   u="<<norm(LbmVec(ux,uy,uz)) <<" RWu="<<norm(LbmVec(rux,ruy,ruz))<<" add"<<doAdd<<" pvel="<<norm(pv) );
                                //if(!bndOk){ mass -= val2fac*size*0.02; } // FORCEDISSOLVE
-                               mass -= val2fac*size*0.0015; // NTEST!
+                               //const LbmFloat val2fac = 1.0; //? TODO N test? /(this->mPartGenProb); // FIXME remove factor!
+                               //mass -= val2fac*size*0.0015; // NTEST!
                                //mass -= val2fac*size*0.001; // NTEST!
+                               mass -= size*0.0020; // NTEST!
 #if LBMDIM==2
                                mpParticles->getLast()->setVel(pv[0],pv[1],0.0);
                                mpParticles->getLast()->setPos(ntlVec3Gfx(srci,srcj,0.5));
 #endif // LBMDIM==2
+               //errMsg("PIT","NEW pit"<<mpParticles->getLast()->getId()<<" pos:"<< mpParticles->getLast()->getPos()<<" status:"<<convertFlags2String(mpParticles->getLast()->getFlags())<<" vel:"<< mpParticles->getLast()->getVel()  );
                                } // multiple parts
                        } // doAdd
                } // */
@@ -2487,17 +2507,19 @@ void LbmFsgrSolver::reinitFlags( int workSet ) {
 
        /* remove empty interface cells that are not allowed to be removed anyway
         * this is important, otherwise the dreaded cell-type-flickering can occur! */
-  for( vector<LbmPoint>::iterator iter=mListEmpty.begin();
-       iter != mListEmpty.end(); iter++ ) {
-    int i=iter->x, j=iter->y, k=iter->z;
+  //for( vector<LbmPoint>::iterator iter=mListEmpty.begin(); iter != mListEmpty.end(); iter++ ) {
+  //int i=iter->x, j=iter->y, k=iter->z;
+       //iter = mListEmpty.erase(iter); iter--; // and continue with next...
+  for(int n=0; n<(int)mListEmpty.size(); n++) {
+    int i=mListEmpty[n].x, j=mListEmpty[n].y, k=mListEmpty[n].z;
                if((RFLAG(workLev,i,j,k, workSet)&(CFInter|CFNoDelete)) == (CFInter|CFNoDelete)) {
-                       // remove entry
-                       if(debugFlagreinit) errMsg("EMPT REMOVED!!!", PRINT_IJK<<" mss"<<QCELL(workLev, i,j,k, workSet, dMass) <<" rho"<< QCELL(workLev, i,j,k, workSet, 0) ); // DEBUG SYMM
-                       iter = mListEmpty.erase(iter); 
-                       iter--; // and continue with next...
-
                        // treat as "new inter"
                        addToNewInterList(i,j,k);
+                       // remove entry
+                       if(debugFlagreinit) errMsg("EMPT REMOVED!!!", PRINT_IJK<<" mss"<<QCELL(workLev, i,j,k, workSet, dMass) <<" rho"<< QCELL(workLev, i,j,k, workSet, 0) ); // DEBUG SYMM
+                       if(n<(int)mListEmpty.size()-1) mListEmpty[n] = mListEmpty[mListEmpty.size()-1]; 
+                       mListEmpty.pop_back();
+                       n--; 
                }
        } 
 
index 589f219..605bc02 100644 (file)
 
 
 
-#if LBM_INCLUDE_TESTSOLVERS!=1
-
 // off for non testing
 #define PRECOLLIDE_MODS(rho,ux,uy,uz) 
 
-#else // LBM_INCLUDE_TESTSOLVERS!=1
-
-#define PRECOLLIDE_GETPOS \
-       LbmVec( \
-                       ((this->mvGeoEnd[0]-this->mvGeoStart[0])/(LbmFloat)this->mSizex) * (LbmFloat)i + this->mvGeoStart[0], \
-                       ((this->mvGeoEnd[1]-this->mvGeoStart[1])/(LbmFloat)this->mSizey) * (LbmFloat)j + this->mvGeoStart[1], \
-                       ((this->mvGeoEnd[2]-this->mvGeoStart[2])/(LbmFloat)this->mSizez) * (LbmFloat)k + this->mvGeoStart[2]  \
-                       )
-
-// debug modifications of collide vars (testing)
-#define PRECOLLIDE_MODS(rho,ux,uy,uz) \
-       if( (this->mTForceStrength>0.) && (RFLAG(0,i,j,k, mLevel[0].setCurr)&CFNoBndFluid) ) { \
-               LbmVec pos = PRECOLLIDE_GETPOS; \
-               LbmVec vel(ux,uy,uz); \
-               mpTest->mControlParts.modifyVelocity(pos,vel); /* */\
-               if((i==16)&&(j==10)){ /*debugMarkCell(0,16,10,0);*/ errMsg("FTDEB"," at "<<PRINT_IJK<<" targ:"<<vel<<",len:"<<norm(vel)<<",  org:"<<ntlVec3Gfx(ux,uy,uz) ); }\
-               ux = vel[0]; uy=vel[1]; uz=vel[2]; \
-               /* test acutal values...? */ \
-       }
-       /* end PRECOLLIDE_MODS */
-
-// debug modifications of collide vars (testing)
-#define _PRECOLLIDE_MODS(rho,ux,uy,uz) \
-       if(this->mTForceStrength>0.) { \
-               LbmVec u(ux,uy,uz); \
-               LbmVec pos = PRECOLLIDE_GETPOS; \
-               LbmVec targv = u; const int lev=0; \
-               mpTest->mControlParts.modifyVelocity(pos,targv); /* */\
-               LbmVec devia = (targv-u); \
-               LbmVec set; \
-               set = u + (devia/this->mOmega) * this->mTForceStrength; /* */\
-               set = u + targv*this->mTForceStrength; /* */\
-               ux = set[0]; uy=set[1]; uz=set[2]; \
-               if(j==10){ errMsg("FTDEB"," at "<<PRINT_IJK<<" targ"<<targv<<","<<norm(targv)<<"  set:"<<set<<","<<norm(set) ); }\
-       }
-       /* end PRECOLLIDE_MODS */
-
-#endif // LBM_INCLUDE_TESTSOLVERS!=1
-
-/*
-
-       \
-       if((0) && (j>=0.01*this->mSizey) && (j<0.9*this->mSizey)) { \
-               if((1) && (i>=0.5*this->mSizex) && (i<0.6*this->mSizex)) { \
-                       LbmFloat len = sqrtf(ux*ux + uy*uy + uz*uz); \
-                       LbmVec u(ux,uy,uz); LbmVec t(1., 0., 0.); \
-                       ux = len*t[0]; uy=len*t[1]; uz=len*t[2]; \
-               } \
-               if((1) && (i>=0.65*this->mSizex) && (i<0.75*this->mSizex)) { \
-                       LbmFloat len = sqrtf(ux*ux + uy*uy + uz*uz); \
-                       LbmVec u(ux,uy,uz); LbmVec t(0., 1., 0.); \
-                       ux = len*t[0]; uy=len*t[1]; uz=len*t[2]; \
-               } \
-       } \
-       \
-       if((0) && (j>=0.0*this->mSizey) && (j<1.0 *this->mSizey)) { \
-               if((1) && (i>=0.6 *this->mSizex) && (i<0.7 *this->mSizex)) { \
-                       LbmFloat len = sqrtf(ux*ux + uy*uy + uz*uz); \
-                       LbmVec u(ux,uy,uz); \
-                       LbmVec t(0., 1., 0.); \
-                       LbmFloat devia = norm(u-t); \
-                       // 0.0001-strong, 0.000001-small, 
-                       t = u - (devia/this->mOmega) * this->mTForceStrength; \
-                       // ux = len*t[0]; uy=len*t[1]; uz=len*t[2]; 
-                       ux = t[0]; uy=t[1]; uz=t[2]; \
-                       //ux= uy= uz= 0.0; 
-                       //errMsg("DDDD"," at "<<PRINT_IJK);  
-               } \
-       } \
-       if(0) { \
-               LbmFloat len = sqrtf(ux*ux + uy*uy + uz*uz); \
-               LbmVec u(ux,uy,uz); \
-               LbmVec p(i/(LbmFloat)this->mSizex, j/(LbmFloat)this->mSizey, k/(LbmFloat)this->mSizez); \
-               LbmVec cp(0.5, 0.5, 0.5); \
-               LbmVec delt = cp - p; \
-               LbmFloat dist = norm(delt); \
-               normalize(delt); \
-               LbmVec tang = cross(delt, LbmVec(0.,0.,1.)); \
-               normalize(tang); \
-               const LbmFloat falloff = 5.0; \
-               LbmVec targv = tang*1.0 * 1.0/(1.0+falloff*dist); \
-               LbmVec devia = (targv-u); \
-               LbmFloat devial = norm(devia); \
-               LbmVec set; \
-               set = u +targv * this->mTForceStrength; \
-               set = u + (devia/this->mOmega) * this->mTForceStrength; \
-               ux = set[0]; uy=set[1]; uz=set[2]; \
-               if(j==10){ errMsg("FTDEB"," at "<<PRINT_IJK<<" dist:"<<delt<<","<<dist<<" targ"<<targv<<","<<norm(targv)<<"  set:"<<set<<","<<norm(set) ); }\
-       } \
-
-*/
        
 /******************************************************************************
  * normal relaxation
@@ -1152,14 +1059,14 @@ inline LbmFloat LbmFsgrSolver::getLesOmega(LbmFloat omega, LbmFloat csmago, LbmF
 }
 
 #define DEBUG_CALCPRINTCELL(str,df) {\
-               LbmFloat rho=df[0], ux=0., uy=0., uz=0.; \
+               LbmFloat prho=df[0], pux=0., puy=0., puz=0.; \
                for(int l=1; l<this->cDfNum; l++) { \
-                       rho += df[l];  \
-                       ux  += (this->dfDvecX[l]*df[l]);  \
-                       uy  += (this->dfDvecY[l]*df[l]);  \
-                       uz  += (this->dfDvecZ[l]*df[l]);  \
+                       prho += df[l];  \
+                       pux  += (this->dfDvecX[l]*df[l]);  \
+                       puy  += (this->dfDvecY[l]*df[l]);  \
+                       puz  += (this->dfDvecZ[l]*df[l]);  \
                } \
-               errMsg("DEBUG_CALCPRINTCELL",">"<<str<<" rho="<<rho<<" vel="<<ntlVec3Gfx(ux,uy,uz) ); \
+               errMsg("DEBUG_CALCPRINTCELL",">"<<str<<" rho="<<prho<<" vel="<<ntlVec3Gfx(pux,puy,puz) ); \
        } /* END DEBUG_CALCPRINTCELL */ 
 
 // "normal" collision
index 7406c55..8cd698e 100644 (file)
@@ -110,39 +110,6 @@ void LbmFsgrSolver::prepareVisualization( void ) {
                *this->mpIso->lbmGetData( i+1 , j+1 ,ZKOFF+ZKD1) += ( val * mIsoWeight[26] ); 
        }
 
-       
-#if LBM_INCLUDE_TESTSOLVERS==1
-       /*if(mUseTestdata) {
-               int border = 1;
-               for(int k=0;k<mLevel[mMaxRefine].lSizez-1;k++)
-                       for(int j=0;j<mLevel[mMaxRefine].lSizey-1;j++) {
-                               for(int l=0; l<=border; l++) {
-                                       *this->mpIso->lbmGetData( l-1,                         j,ZKOFF) = *this->mpIso->lbmGetData( border+1,                         j,ZKOFF);
-                                       *this->mpIso->lbmGetData( mLevel[mMaxRefine].lSizex-l, j,ZKOFF) = *this->mpIso->lbmGetData( mLevel[mMaxRefine].lSizex-border-1, j,ZKOFF);
-                               }
-                       }
-
-               for(int k=0;k<mLevel[mMaxRefine].lSizez-1;k++)
-                       for(int i=-1;i<mLevel[mMaxRefine].lSizex+1;i++) {      
-                               for(int l=0; l<=border; l++) {
-                                       *this->mpIso->lbmGetData( i, l-1,                         ZKOFF) = *this->mpIso->lbmGetData( i, border+1,                         ZKOFF);
-                                       *this->mpIso->lbmGetData( i, mLevel[mMaxRefine].lSizey-l, ZKOFF) = *this->mpIso->lbmGetData( i, mLevel[mMaxRefine].lSizey-border-1, ZKOFF);
-                               }
-                       }
-
-               if(LBMDIM == 3) {
-                       // only for 3D
-                       for(int j=-1;j<mLevel[mMaxRefine].lSizey+1;j++)
-                               for(int i=-1;i<mLevel[mMaxRefine].lSizex+1;i++) {      
-                                       for(int l=0; l<=border; l++) {
-                                               *this->mpIso->lbmGetData( i,j,l-1                        ) = *this->mpIso->lbmGetData( i,j, border+1                         );
-                                               *this->mpIso->lbmGetData( i,j,mLevel[mMaxRefine].lSizez-l) = *this->mpIso->lbmGetData( i,j,mLevel[mMaxRefine].lSizez-1-border);
-                                       }
-                               }
-               }
-       } // testdata */
-#endif // LBM_INCLUDE_TESTSOLVERS==1
-       // */
 
        // update preview, remove 2d?
        if((this->mOutputSurfacePreview)&&(LBMDIM==3)) {
@@ -175,7 +142,7 @@ void LbmFsgrSolver::prepareVisualization( void ) {
                                *mpPreviewSurface->lbmGetData(i,j,pvsz-1) = *this->mpIso->lbmGetData( (int)(i*scalex), (int)(j*scaley) , this->mSizez-1);
                        }
 
-               if(mUseTestdata) {
+               if(mFarFieldSize>=2.) {
                        // also remove preview border
                        for(int k= 0; k< (pvsz-1); ++k) {
                                for(int j=0;j< pvsy;j++) {
@@ -267,11 +234,11 @@ vector<ntlGeometryObject*> LbmFsgrSolver::getDebugObjects() {
  *****************************************************************************/
 
 /*! init particle positions */
-int LbmFsgrSolver::initParticles(ParticleTracer *partt) { 
+int LbmFsgrSolver::initParticles() { 
   int workSet = mLevel[mMaxRefine].setCurr;
   int tries = 0;
   int num = 0;
-       mpParticles=partt;
+       ParticleTracer *partt = mpParticles;
 
   partt->setStart( this->mvGeoStart + ntlVec3Gfx(mLevel[mMaxRefine].nodeSize*0.5) );
   partt->setEnd  ( this->mvGeoEnd   + ntlVec3Gfx(mLevel[mMaxRefine].nodeSize*0.5) );
@@ -279,11 +246,14 @@ int LbmFsgrSolver::initParticles(ParticleTracer *partt) {
   partt->setSimStart( ntlVec3Gfx(0.0) );
   partt->setSimEnd  ( ntlVec3Gfx(this->mSizex,   this->mSizey,   getForZMaxBnd(mMaxRefine)) );
   
-  while( (num<partt->getNumParticles()) && (tries<100*partt->getNumParticles()) ) {
+  while( (num<partt->getNumInitialParticles()) && (tries<100*partt->getNumInitialParticles()) ) {
     LbmFloat x,y,z;
-    x = 0.0+(( (LbmFloat)(this->mSizex-1) )     * (rand()/(RAND_MAX+1.0)) );
-    y = 0.0+(( (LbmFloat)(this->mSizey-1) )     * (rand()/(RAND_MAX+1.0)) );
-    z = 0.0+(( (LbmFloat) getForZMax1(mMaxRefine) )* (rand()/(RAND_MAX+1.0)) );
+    //x = 0.0+(( (LbmFloat)(this->mSizex-1) )     * (rand()/(RAND_MAX+1.0)) );
+    //y = 0.0+(( (LbmFloat)(this->mSizey-1) )     * (rand()/(RAND_MAX+1.0)) );
+    //z = 0.0+(( (LbmFloat) getForZMax1(mMaxRefine) )* (rand()/(RAND_MAX+1.0)) );
+    x = 1.0+(( (LbmFloat)(this->mSizex-3.) )          * (rand()/(RAND_MAX+1.0)) );
+    y = 1.0+(( (LbmFloat)(this->mSizey-3.) )          * (rand()/(RAND_MAX+1.0)) );
+    z = 1.0+(( (LbmFloat) getForZMax1(mMaxRefine)-2. )* (rand()/(RAND_MAX+1.0)) );
     int i = (int)(x+0.5);
     int j = (int)(y+0.5);
     int k = (int)(z+0.5);
@@ -292,12 +262,19 @@ int LbmFsgrSolver::initParticles(ParticleTracer *partt) {
       z = 0.5; // place in the middle of domain
     }
 
-    if( TESTFLAG( RFLAG(mMaxRefine, i,j,k, workSet), CFFluid ) ||
-        TESTFLAG( RFLAG(mMaxRefine, i,j,k, workSet), CFInter ) ) { // only fluid cells?
+    //if( TESTFLAG( RFLAG(mMaxRefine, i,j,k, workSet), CFFluid ) ||
+        //TESTFLAG( RFLAG(mMaxRefine, i,j,k, workSet), CFInter ) ) { // only fluid cells?
+    if( TESTFLAG( RFLAG(mMaxRefine, i,j,k, workSet), CFFluid ) 
+        //&& TESTFLAG( RFLAG(mMaxRefine, i,j,k, workSet), CFNoNbFluid ) 
+                               ) { // inner fluid only 
+                       bool cellOk = true;
+                       FORDF1 { if(!(RFLAG_NB(mMaxRefine,i,j,k,workSet, l) & CFFluid)) cellOk = false; }
+                       if(!cellOk) continue;
       // in fluid...
       partt->addParticle(x,y,z);
                        partt->getLast()->setStatus(PART_IN);
-                       partt->getLast()->setType(PART_BUBBLE);
+                       partt->getLast()->setType(PART_TRACER);
+                       partt->getLast()->setSize(1.0);
       num++;
     }
     tries++;
@@ -337,6 +314,7 @@ int LbmFsgrSolver::initParticles(ParticleTracer *partt) {
                                if(partDebug) errMsg("PARTTT","SET "<<PRINT_VEC(x,y,z)<<" p"<<partt->getLast()->getPos() <<" s"<<partt->getLast()->getSize() );
                        }
                }
+               // place floats on rectangular region FLOAT_JITTER_BND
                if(mpTest->mDebugvalue2==-10.0){ 
                        const int lev = mMaxRefine;
                        const int sx = mLevel[lev].lSizex;
@@ -348,8 +326,9 @@ int LbmFsgrSolver::initParticles(ParticleTracer *partt) {
                                        LbmFloat x,y,z;
                                        x = 0.0+(LbmFloat)(i);
                                        y = 0.0+(LbmFloat)(j);
-                                       //z = 0.5+(LbmFloat)(k);
-                                       z = mLevel[lev].lSizez/20.0*8.0 - 1.0;
+                                       //z = mLevel[lev].lSizez/10.0*2.5 - 1.0;
+                                       z = mLevel[lev].lSizez/20.0*7.5 - 1.0;
+                                       //z = mLevel[lev].lSizez/20.0*4.5 - 1.0;
                                        partt->addParticle(x,y,z);
                                        //if( (i>0)&&(i<sx) && (j>0)&&(j<sy) ) { partt->getLast()->setStatus(PART_IN); } else { partt->getLast()->setStatus(PART_OUT); }
                                        partt->getLast()->setStatus(PART_IN);
@@ -363,17 +342,21 @@ int LbmFsgrSolver::initParticles(ParticleTracer *partt) {
 #endif // LBM_INCLUDE_TESTSOLVERS
 
        
-  debMsgStd("LbmFsgrSolver::initParticles",DM_MSG,"Added "<<num<<" particles, genProb:"<<this->mPartGenProb, 10);
+  debMsgStd("LbmFsgrSolver::initParticles",DM_MSG,"Added "<<num<<" particles, genProb:"<<this->mPartGenProb<<", tries:"<<tries, 10);
   if(num != partt->getNumParticles()) return 1;
 
        return 0;
 }
 
-void LbmFsgrSolver::advanceParticles(ParticleTracer *partt) { 
+#define P_CHANGETYPE(p, newtype) \
+               p->setLifeTime(0); \
+    /* errMsg("PIT","U pit"<<(p)->getId()<<" pos:"<< (p)->getPos()<<" status:"<<convertFlags2String((p)->getFlags())<<" to "<< (newtype) ); */ \
+               p->setType(newtype); 
+
+void LbmFsgrSolver::advanceParticles() { 
   int workSet = mLevel[mMaxRefine].setCurr;
        LbmFloat vx=0.0,vy=0.0,vz=0.0;
        LbmFloat rho, df[27]; //feq[27];
-       if(mpParticles!=partt) { errMsg("LbmFsgrSolver::advanceParticles","Invalid ParticleTracer..."); }
 #define DEL_PART { \
        /*errMsg("PIT","DEL AT "<< __LINE__<<" type:"<<p->getType()<<"  ");  */ \
        p->setActive( false ); \
@@ -394,6 +377,7 @@ void LbmFsgrSolver::advanceParticles(ParticleTracer *partt) {
        const LbmFloat v1 = 9.0; // v max
        const LbmFloat v2 = 2.0; // v min
        const LbmVec rwgrav = vec2L( this->mpParam->getGravity(mSimulationTime) );
+       const bool useff = (mFarFieldSize>2.); // if(mpTest->mDebugvalue1>0.0){ 
 
        // TODO use timestep size
        //bool isIn,isOut,isInZ;
@@ -402,11 +386,11 @@ void LbmFsgrSolver::advanceParticles(ParticleTracer *partt) {
        //const int cutval = 0; // TODO FIXME add half border!
        const int cutval = mCutoff; // use full border!?
        int actCnt=0;
-       if(this->mStepCnt%50==49) { partt->cleanup(); }
-  for(vector<ParticleObject>::iterator pit= partt->getParticlesBegin();
-      pit!= partt->getParticlesEnd(); pit++) {
-    //errMsg("PIT"," pit "<< (*pit).getPos()<<" status:"<<convertFlags2String((*pit).getFlags())<<" vel:"<< (*pit).getVel()  );
-    //errMsg("PIT"," pit pos:"<< (*pit).getPos()<<" vel:"<< (*pit).getVel()<<" status:"<<convertFlags2String((*pit).getFlags()) <<" " <<partt->getStart()<<" "<<partt->getEnd() );
+       if(this->mStepCnt%50==49) { mpParticles->cleanup(); }
+  for(vector<ParticleObject>::iterator pit= mpParticles->getParticlesBegin();
+      pit!= mpParticles->getParticlesEnd(); pit++) {
+    //errMsg("PIT"," pit"<<(*pit)->getId()<<" pos:"<< (*pit).getPos()<<" status:"<<convertFlags2String((*pit).getFlags())<<" vel:"<< (*pit).getVel()  );
+    //errMsg("PIT"," pit pos:"<< (*pit).getPos()<<" vel:"<< (*pit).getVel()<<" status:"<<convertFlags2String((*pit).getFlags()) <<" " <<mpParticles->getStart()<<" "<<mpParticles->getEnd() );
                //int flag = (*pit).getFlags();
     if( (*pit).getActive()==false ) continue;
     int i,j,k;
@@ -418,25 +402,23 @@ void LbmFsgrSolver::advanceParticles(ParticleTracer *partt) {
                i= (int)(pos[0]+0.5);
                j= (int)(pos[1]+0.5);
                k= (int)(pos[2]+0.5);
-               if(LBMDIM==2) {
-                       k = 0;
-               }
+               if(LBMDIM==2) { k = 0; }
 
                // only testdata handling, all for sws
 #if LBM_INCLUDE_TESTSOLVERS==1
-               if(mUseTestdata) { 
-                       if(mpTest->mDebugvalue1>0.0){ 
-                               p->setStatus(PART_OUT);
-                               mpTest->handleParticle(p, i,j,k); continue;
-               } }
+               if(useff) {
+                       p->setStatus(PART_OUT);
+                       mpTest->handleParticle(p, i,j,k); continue;
+               } 
 #endif // LBM_INCLUDE_TESTSOLVERS==1
 
+               // FIXME , add k tests again, remove per type ones...
                if(p->getStatus()&PART_IN) { // IN
                        if( (i<cutval)||(i>this->mSizex-1-cutval)||
                                        (j<cutval)||(j>this->mSizey-1-cutval)
                                        //||(k<cutval)||(k>this->mSizez-1-cutval) 
                                        ) {
-                               if(!mUseTestdata) { DEL_PART;
+                               if(!useff) { DEL_PART;
                                } else { 
                                        p->setStatus(PART_OUT); 
                                        /* del? */ //if((rand()/(RAND_MAX+1.0))<0.5) DEL_PART;
@@ -454,7 +436,8 @@ void LbmFsgrSolver::advanceParticles(ParticleTracer *partt) {
                }
                //p->setStatus(PART_OUT);// DEBUG always out!
 
-               if(p->getType()==PART_BUBBLE) {
+               if( (p->getType()==PART_BUBBLE) ||
+                   (p->getType()==PART_TRACER) ) {
 
                        // no interpol
                        rho = vx = vy = vz = 0.0;
@@ -481,14 +464,14 @@ void LbmFsgrSolver::advanceParticles(ParticleTracer *partt) {
                                        } else { // OUT
                                                // out of bounds, deactivate...
                                                // FIXME make fsgr treatment
-                                               p->setType( PART_FLOAT ); continue;
+                                               if(p->getType()==PART_BUBBLE) { P_CHANGETYPE(p, PART_FLOAT ); continue; }
                                        }
                                } else {
                                        // below 3d region, just rise
                                }
                        } else { // OUT
 #if LBM_INCLUDE_TESTSOLVERS==1
-                               if(mUseTestdata) { mpTest->handleParticle(p, i,j,k); }
+                               if(useff) { mpTest->handleParticle(p, i,j,k); }
                                else DEL_PART;
 #else // LBM_INCLUDE_TESTSOLVERS==1
                                DEL_PART;
@@ -497,7 +480,7 @@ void LbmFsgrSolver::advanceParticles(ParticleTracer *partt) {
                        }
 
                        ntlVec3Gfx v = p->getVel(); // dampen...
-                       if(mUseTestdata) {
+                       if( (useff)&& (p->getType()==PART_BUBBLE) ) {
                                // test rise
                                //O vz = p->getVel()[2]-0.5*mLevel[mMaxRefine].gravity[2];
 
@@ -533,11 +516,35 @@ void LbmFsgrSolver::advanceParticles(ParticleTracer *partt) {
                                LbmFloat w = 0.99;
                                vz = (1.0-w)*vz + w*(p->getVel()[2]-0.5*(p->getSize()/5.0)*mLevel[mMaxRefine].gravity[2]);
                                v = ntlVec3Gfx(vx,vy,vz)+vec2G(fd2);
+                       } else if(p->getType()==PART_TRACER) {
+                               v = ntlVec3Gfx(vx,vy,vz);
+                               if( RFLAG(mMaxRefine, i,j,k, workSet)&(CFFluid) ) {
+                                       // ok
+                               } else {
+                                       const int lev = mMaxRefine;
+                                       LbmFloat nx,ny,nz, nv1,nv2;
+                       //mynbfac = QCELL_NB(lev, i,j,k,SRCS(lev),l, dFlux) / QCELL(lev, i,j,k,SRCS(lev), dFlux);
+                                       if(RFLAG_NB(lev,i,j,k,workSet, dE) &(CFFluid|CFInter)){ nv1 = QCELL_NB(lev,i,j,k,workSet,dE,dFfrac); } else nv1 = 0.0;
+                                       if(RFLAG_NB(lev,i,j,k,workSet, dW) &(CFFluid|CFInter)){ nv2 = QCELL_NB(lev,i,j,k,workSet,dW,dFfrac); } else nv2 = 0.0;
+                                       nx = 0.5* (nv2-nv1);
+                                       if(RFLAG_NB(lev,i,j,k,workSet, dN) &(CFFluid|CFInter)){ nv1 = QCELL_NB(lev,i,j,k,workSet,dN,dFfrac); } else nv1 = 0.0;
+                                       if(RFLAG_NB(lev,i,j,k,workSet, dS) &(CFFluid|CFInter)){ nv2 = QCELL_NB(lev,i,j,k,workSet,dS,dFfrac); } else nv2 = 0.0;
+                                       ny = 0.5* (nv2-nv1);
+#if LBMDIM==3
+                                       if(RFLAG_NB(lev,i,j,k,workSet, dT) &(CFFluid|CFInter)){ nv1 = QCELL_NB(lev,i,j,k,workSet,dT,dFfrac); } else nv1 = 0.0;
+                                       if(RFLAG_NB(lev,i,j,k,workSet, dB) &(CFFluid|CFInter)){ nv2 = QCELL_NB(lev,i,j,k,workSet,dB,dFfrac); } else nv2 = 0.0;
+                                       nz = 0.5* (nv2-nv1);
+#else //LBMDIM==3
+                                       nz = 0.0;
+#endif //LBMDIM==3
+
+                                       v = (ntlVec3Gfx(nx,ny,nz)) * -0.025;
+                               }
                        }
+
                        p->setVel( v );
-                       //p->setVel( ntlVec3Gfx(vx,vy,vz) );
                        p->advanceVel();
-                 // fluid particle
+                       //errMsg("PPPP"," pos"<<p->getPos()<<" "<<p->getVel() );
                } 
 
                // drop handling
@@ -583,13 +590,17 @@ void LbmFsgrSolver::advanceParticles(ParticleTracer *partt) {
                                if(k<cutval) { DEL_PART; continue; }
                                if(k<=this->mSizez-1-cutval){ 
                                        //if( RFLAG(mMaxRefine, i,j,k, workSet)& (CFEmpty|CFInter)) {
-                                       if( RFLAG(mMaxRefine, i,j,k, workSet)& (CFEmpty)) {
+                                       if( RFLAG(mMaxRefine, i,j,k, workSet)& (CFEmpty|CFInter|CFBnd)) {
                                                // still ok
-                                       } else if( RFLAG(mMaxRefine, i,j,k, workSet) & (CFFluid|CFInter) ){
+                                       // shipt3 } else if( RFLAG(mMaxRefine, i,j,k, workSet) & (CFFluid|CFInter) ){
+                                       } else if( RFLAG(mMaxRefine, i,j,k, workSet) & (CFFluid) ){
                                                // FIXME make fsgr treatment
-                                               if(p->getLifeTime()>50) { 
-                                                       p->setType( PART_FLOAT ); continue; 
-                                               } else DEL_PART;
+                                               //if(p->getLifeTime()>50) { 
+                                               P_CHANGETYPE(p, PART_FLOAT ); continue; 
+                                               // jitter in cell to prevent stacking when hitting a steep surface
+                                               LbmVec pos = p->getPos(); pos[0] += (rand()/(RAND_MAX+1.0))-0.5;
+                                               pos[1] += (rand()/(RAND_MAX+1.0))-0.5; p->setPos(pos);
+                                               //} else DEL_PART;
                                        } else {
                                                DEL_PART;
                                                this->mNumParticlesLost++;
@@ -597,7 +608,7 @@ void LbmFsgrSolver::advanceParticles(ParticleTracer *partt) {
                                }
                        } else { // OUT
 #if LBM_INCLUDE_TESTSOLVERS==1
-                               if(mUseTestdata) { mpTest->handleParticle(p, i,j,k); }
+                               if(useff) { mpTest->handleParticle(p, i,j,k); }
                                else{ DEL_PART; }
 #else // LBM_INCLUDE_TESTSOLVERS==1
                                { DEL_PART; }
@@ -619,9 +630,10 @@ void LbmFsgrSolver::advanceParticles(ParticleTracer *partt) {
                                } else if( TESTFLAG( RFLAG(mMaxRefine, i,j,k, workSet), CFFluid ) ) {
                        //errMsg("PIT","NEWBUB pit "<< (*pit).getPos()<<" status:"<<convertFlags2String((*pit).getFlags())  );
 
-                                       //p->setType( PART_BUBBLE ); continue;
+                                       //P_CHANGETYPE(p, PART_BUBBLE ); continue;
                                        // currently bubbles off! DEBUG!
-                                       DEL_PART; // DEBUG bubbles off for siggraph
+                                       //DEL_PART; // DEBUG bubbles off for siggraph
+                                       P_CHANGETYPE(p, PART_FLOAT ); continue;
                                } else if( TESTFLAG( RFLAG(mMaxRefine, i,j,k, workSet), CFEmpty ) ) {
                        //errMsg("PIT","NEWDROP pit "<< (*pit).getPos()<<" status:"<<convertFlags2String((*pit).getFlags())  );
                                        //? if(p->getLifeTime()>50) {
@@ -631,13 +643,14 @@ void LbmFsgrSolver::advanceParticles(ParticleTracer *partt) {
                                                        //(j<=cutval)||(j>=this->mSizey-1-cutval)) {
                                        //} else 
                                        //if(p->getLifeTime()>10) {
-                                       p->setType( PART_DROP ); continue;
+                                       P_CHANGETYPE(p, PART_DROP ); continue;
                                        //} else DEL_PART;                                              
                                        
                                }
                        } else { // OUT
                                // undecided particle outside... remove?
                                DEL_PART; 
+                               //P_CHANGETYPE(p, PART_FLOAT ); continue;
                        }
                }
 
@@ -650,11 +663,17 @@ void LbmFsgrSolver::advanceParticles(ParticleTracer *partt) {
 #if LBM_INCLUDE_TESTSOLVERS==1
                        // vanishing
                        prob = (rand()/(RAND_MAX+1.0));
-                       if((mUseTestdata)&&(k>mpTest->mFluidHeight)) {
-                               LbmFloat fac = (LbmFloat)(k-mpTest->mFluidHeight)/(LbmFloat)(10*(mLevel[mMaxRefine].lSizez-mpTest->mFluidHeight));
-                               prob /= fac; //  TODO test? errMsg("T","T "<<prob<<" "<<fac);
+                       // increse del prob. up to max height by given factor
+                       const int fhCheckh = (int)(mpTest->mFluidHeight*1.25);
+                       const LbmFloat fhProbh = 25.;
+                       if((useff)&&(k>fhCheckh)) {
+                               LbmFloat fac = (LbmFloat)(k-fhCheckh)/(LbmFloat)(fhProbh*(mLevel[mMaxRefine].lSizez-fhCheckh));
+                               prob /= fac; 
                        }
                        if(prob<mLevel[mMaxRefine].timestep*0.1) DEL_PART;
+                       // forced vanishing
+                       //? if(k>this->mSizez*3/4) {    if(prob<3.0*mLevel[mMaxRefine].timestep*0.1) DEL_PART;}
+
 #else // LBM_INCLUDE_TESTSOLVERS==1
 #endif // LBM_INCLUDE_TESTSOLVERS==1
 
@@ -663,13 +682,15 @@ void LbmFsgrSolver::advanceParticles(ParticleTracer *partt) {
                                //ntlVec3Gfx v = getVelocityAt(i,j,k);
                                rho = vx = vy = vz = 0.0;
 
-                               //const int DEPTH_AVG=11; // TODO how much!?
-                               const int DEPTH_AVG=7; // TODO how much!?
+                               // define from particletracer.h
+#if MOVE_FLOATS==1
+                               const int DEPTH_AVG=3; // only average interface vels
                                int ccnt=0;
-                               for(int kk=1;kk<DEPTH_AVG;kk+=2) {
+                               for(int kk=0;kk<DEPTH_AVG;kk+=1) {
                                //for(int kk=1;kk<DEPTH_AVG;kk+=1) {
                                        if((k-kk)<1) continue;
-                                       if(RFLAG(mMaxRefine, i,j,k, workSet)&(CFFluid|CFInter)) {} else continue;
+                                       //if(RFLAG(mMaxRefine, i,j,k, workSet)&(CFFluid|CFInter)) {} else continue;
+                                       if(RFLAG(mMaxRefine, i,j,k, workSet)&(CFInter)) {} else continue;
                                        ccnt++;
                                        FORDF0{
                                                LbmFloat cdf = QCELL(mMaxRefine, i,j,k-kk, workSet, l);
@@ -681,32 +702,44 @@ void LbmFsgrSolver::advanceParticles(ParticleTracer *partt) {
                                        }
                                }
                                if(ccnt) {
-                               vx /=(LbmFloat)(ccnt * 1.0); // half xy speed! value2
-                               vy /=(LbmFloat)(ccnt * 1.0);
+                               // use halved surface velocity (todo, use omega instead)
+                               vx /=(LbmFloat)(ccnt * 2.0); // half xy speed! value2
+                               vy /=(LbmFloat)(ccnt * 2.0);
                                vz /=(LbmFloat)(ccnt); }
-                               // forced vanishing
-                               if(k>this->mSizez*3/4) {        if(prob<3.0*mLevel[mMaxRefine].timestep*0.1) DEL_PART;}
+#else // MOVE_FLOATS==1
+                               vx=vy=0.; //p->setVel(ntlVec3Gfx(0.) ); // static_float
+#endif // MOVE_FLOATS==1
+                               vx += (rand()/(RAND_MAX+1.0))*FLOAT_JITTER-(FLOAT_JITTER*0.5);
+                               vy += (rand()/(RAND_MAX+1.0))*FLOAT_JITTER-(FLOAT_JITTER*0.5);
 
+                               bool delfloat = false;
                                if( TESTFLAG( RFLAG(mMaxRefine, i,j,k, workSet), CFFluid ) ) {
+                                       // in fluid cell
                                        if((1) && (k<this->mSizez-3) && 
                                                        (
                                                          TESTFLAG( RFLAG(mMaxRefine, i,j,k+1, workSet), CFInter ) ||
-                                                         TESTFLAG( RFLAG(mMaxRefine, i,j,k+2, workSet), CFInter ) )
-                                                        ) {
+                                                         TESTFLAG( RFLAG(mMaxRefine, i,j,k+2, workSet), CFInter ) 
+                                                       ) ) {
                                                vz = p->getVel()[2]-0.5*mLevel[mMaxRefine].gravity[2];
                                                if(vz<0.0) vz=0.0;
-                                       } else DEL_PART;
+                                       } else delfloat=true;
+                                       // keep below obstacles
+                                       if((delfloat) && (TESTFLAG( RFLAG(mMaxRefine, i,j,k+1, workSet), CFBnd )) ) {
+                                               delfloat=false; vz=0.0;
+                                       }
                                } else if( TESTFLAG( RFLAG(mMaxRefine, i,j,k, workSet), CFInter ) ) {
                                        // keep in interface , one grid cell offset is added in part. gen
-                               } 
-                               // check if above inter, remove otherwise
-                               else if((1) && (k>2) && (
-                                                       TESTFLAG( RFLAG(mMaxRefine, i,j,k-1, workSet), CFInter ) ||
-                                                       TESTFLAG( RFLAG(mMaxRefine, i,j,k-2, workSet), CFInter ) )
-                                                       ) {
-                                       vz = p->getVel()[2]+0.5*mLevel[mMaxRefine].gravity[2];
-                                       if(vz>0.0) vz=0.0;
-                               } else DEL_PART; // */
+                               } else { //if( TESTFLAG( RFLAG(mMaxRefine, i,j,k, workSet), CFEmpty ) ) { // shipt?, was all else before
+                                       // check if above inter, remove otherwise
+                                       if((1) && (k>2) && (
+                                                               TESTFLAG( RFLAG(mMaxRefine, i,j,k-1, workSet), CFInter ) ||
+                                                               TESTFLAG( RFLAG(mMaxRefine, i,j,k-2, workSet), CFInter ) 
+                                                               ) ) {
+                                               vz = p->getVel()[2]+0.5*mLevel[mMaxRefine].gravity[2];
+                                               if(vz>0.0) vz=0.0;
+                                       } else delfloat=true; // */
+                               }
+                               if(delfloat) DEL_PART;
                                /*
                                // move down from empty
                                else if( TESTFLAG( RFLAG(mMaxRefine, i,j,k, workSet), CFEmpty ) ) {
@@ -715,20 +748,38 @@ void LbmFsgrSolver::advanceParticles(ParticleTracer *partt) {
                                        //DEL_PART; // ????
                                } else  {        DEL_PART; } // */
                                //vz = 0.0; // DEBUG
-                               ntlVec3Gfx v(vx,vy,vz);
-                               p->setVel( vec2G(v) ); //?
+
+                               p->setVel( vec2G( ntlVec3Gfx(vx,vy,vz) ) ); //?
                                //p->setVel( vec2G(v)*0.75 + p->getVel()*0.25 ); //?
                                p->advanceVel();
                                //errMsg("PIT","IN pit "<< (*pit).getPos()<<" status:"<<convertFlags2String((*pit).getFlags())  );
                        } else {
 #if LBM_INCLUDE_TESTSOLVERS==1
-                               if(mUseTestdata) { mpTest->handleParticle(p, i,j,k); }
+                               if(useff) { mpTest->handleParticle(p, i,j,k); }
                                else DEL_PART; 
 #else // LBM_INCLUDE_TESTSOLVERS==1
                                DEL_PART; 
 #endif // LBM_INCLUDE_TESTSOLVERS==1
                                //errMsg("PIT","OUT pit "<< (*pit).getPos()<<" status:"<<convertFlags2String((*pit).getFlags())  );
                        }
+                               
+                       // additional bnd jitter
+                       if((1) && (useff) && (p->getLifeTime()<3)) {
+                               // use half butoff border 1/8
+                               int maxdw = (int)(mLevel[mMaxRefine].lSizex*0.125*0.5);
+                               if(maxdw<3) maxdw=3;
+#define FLOAT_JITTER_BND (FLOAT_JITTER*2.0)
+#define FLOAT_JITTBNDRAND(x) ((rand()/(RAND_MAX+1.0))*FLOAT_JITTER_BND*(1.-(x/(LbmFloat)maxdw))-(FLOAT_JITTER_BND*(1.-(x)/(LbmFloat)maxdw)*0.5)) 
+                               //if(ABS(i-(               cutval))<maxdw) { p->advance(  (rand()/(RAND_MAX+1.0))*FLOAT_JITTER_BND-(FLOAT_JITTER_BND*0.5)   ,0.,0.); }
+                               //if(ABS(i-(this->mSizex-1-cutval))<maxdw) { p->advance(  (rand()/(RAND_MAX+1.0))*FLOAT_JITTER_BND-(FLOAT_JITTER_BND*0.5)   ,0.,0.); }
+                               if((j>=0)&&(j<=this->mSizey-1)) {
+                                       if(ABS(i-(               cutval))<maxdw) { p->advance(  FLOAT_JITTBNDRAND( ABS(i-(               cutval))), 0.,0.); }
+                                       if(ABS(i-(this->mSizex-1-cutval))<maxdw) { p->advance(  FLOAT_JITTBNDRAND( ABS(i-(this->mSizex-1-cutval))), 0.,0.); }
+                               }
+//#undef FLOAT_JITTER_BND
+#undef FLOAT_JITTBNDRAND
+                               //if( (i<cutval)||(i>this->mSizex-1-cutval)|| //(j<cutval)||(j>this->mSizey-1-cutval)
+                       }
                } 
                
                // unknown particle type        
@@ -737,12 +788,13 @@ void LbmFsgrSolver::advanceParticles(ParticleTracer *partt) {
                }
   }
        myTime_t parttend = getTime(); 
-       debMsgStd("LbmFsgrSolver::advanceParticles",DM_MSG,"Time for particle update:"<< getTimeString(parttend-parttstart)<<" "<<partt->getNumParticles() , 10 );
+       debMsgStd("LbmFsgrSolver::advanceParticles",DM_MSG,"Time for particle update:"<< getTimeString(parttend-parttstart)<<" "<<mpParticles->getNumParticles() , 10 );
 }
 
-void LbmFsgrSolver::notifySolverOfDump(int frameNr,char *frameNrStr,string outfilename) {
+void LbmFsgrSolver::notifySolverOfDump(int dumptype, int frameNr,char *frameNrStr,string outfilename) {
        // debug - raw dump of ffrac values
        if(getenv("ELBEEM_RAWDEBUGDUMP")) {
+               int workSet = mLevel[mMaxRefine].setCurr;
                std::ostringstream name;
                //name <<"fill_" << this->mStepCnt <<".dump";
                name << outfilename<< frameNrStr <<".dump";
@@ -752,7 +804,9 @@ void LbmFsgrSolver::notifySolverOfDump(int frameNr,char *frameNrStr,string outfi
                        for(int k= getForZMinBnd(); k< getForZMaxBnd(mMaxRefine); ++k)  {
                                for(int j=0;j<mLevel[mMaxRefine].lSizey-0;j++)  {
                                        for(int i=0;i<mLevel[mMaxRefine].lSizex-0;i++) {
-                                               float val = QCELL(mMaxRefine,i,j,k, mLevel[mMaxRefine].setCurr,dFfrac);
+                                               float val = 0.;
+                                               if(RFLAG(mMaxRefine, i,j,k, workSet) & CFInter) val = QCELL(mMaxRefine,i,j,k, mLevel[mMaxRefine].setCurr,dFfrac);
+                                               if(RFLAG(mMaxRefine, i,j,k, workSet) & CFFluid) val = 1.;
                                                //fwrite( &val, sizeof(val), 1, file); // binary
                                                fprintf(file, "%f ",val); // text
                                                //errMsg("W", PRINT_IJK<<" val:"<<val);
index 34a2150..24a1889 100644 (file)
@@ -2392,7 +2392,7 @@ void readVelgz(char *filename, Object *srcob)
                        vverts[i].co[j] = 0.; 
                } 
        } 
-       if(srcob->fluidsimSettings->typeFlags&OB_FSDOMAIN_NOVECGEN) return;
+       if(srcob->fluidsimSettings->domainNovecgen>0) return;
 
        if(len<7) { 
                //printf("readVelgz Eror: invalid filename '%s'\n",filename); // DEBUG
@@ -2516,6 +2516,18 @@ void loadFluidsimMesh(Object *srcob, int useRenderParams)
                mesh = readBobjgz(targetFile, srcob->fluidsimSettings->orgMesh, bbSize,bbSize );
        }
        if(!mesh) {
+               // switch, abort background rendering when fluidsim mesh is missing
+               const char *strEnvName2 = "BLENDER_ELBEEMBOBJABORT"; // from blendercall.cpp
+               if(G.background==1) {
+                       if(getenv(strEnvName2)) {
+                               int elevel = atoi(getenv(strEnvName2));
+                               if(elevel>0) {
+                                       printf("Env. var %s set, fluid sim mesh not found, aborting render...\n",strEnvName2);
+                                       exit(1);
+                               }
+                       }
+               }
+               
                // display org. object upon failure
                srcob->data = srcob->fluidsimSettings->orgMesh;
                return;
index 2a95bd8..e0c7509 100644 (file)
@@ -1578,6 +1578,11 @@ static pMatrixCache *cache_object_matrices(Object *ob, int start, int end)
        return mcache;
 }
 
+/* for fluidsim win32 debug messages */
+#if defined(WIN32) && (!(defined snprintf))
+#define snprintf _snprintf
+#endif
+
 /* main particle building function 
    one day particles should become dynamic (realtime) with the current method as a 'bake' (ton) */
 void build_particle_system(Object *ob)
@@ -1599,7 +1604,7 @@ void build_particle_system(Object *ob)
        float *volengths= NULL, *folengths= NULL;
        int deform=0, a, totpart, paf_sta, paf_end;
        int waitcursor_set= 0, totvert, totface, curface, curvert;
-       int readMask =0, activeParts;
+       int readMask, activeParts, fileParts;
        
        /* return conditions */
        if(ob->type!=OB_MESH) return;
@@ -1617,6 +1622,7 @@ void build_particle_system(Object *ob)
                char *suffix  = "fluidsurface_particles_#";
                char *suffix2 = ".gz";
                char filename[256];
+               char debugStrBuffer[256];
                int  curFrame = G.scene->r.cfra -1; // warning - sync with derived mesh fsmesh loading
                int  j, numFileParts;
                gzFile gzf;
@@ -1635,8 +1641,7 @@ void build_particle_system(Object *ob)
 
                gzf = gzopen(filename, "rb");
                if (!gzf) {
-                       //char debugStrBuffer[256];
-                       //define win32... snprintf(debugStrBuffer,256,"readFsPartData::error - Unable to open file for reading '%s'\n", filename);
+                       snprintf(debugStrBuffer,256,"readFsPartData::error - Unable to open file for reading '%s' \n", filename); 
                        //elbeemDebugOut(debugStrBuffer);
                        paf->totpart = 1;
                        return;
@@ -1648,28 +1653,24 @@ void build_particle_system(Object *ob)
                paf->totpart= totpart;
                paf->totkey= 1;
                /* initialize particles */
-               new_particle(paf);// ?
-               ftime = 0.0;
-               dtime= 0.0f;
+               new_particle(paf); 
+               ftime = 0.0; // unused...
 
                // set up reading mask
-               //for(j=1; j<=4; j++ ){ if(ob->fluidsimSettings->guiDisplayMode&j) readMask |= (1<<j); }
-               readMask = ob->fluidsimSettings->guiDisplayMode;
+               readMask = ob->fluidsimSettings->typeFlags;
                activeParts=0;
-               // FIXME only allocate needed ones?
+               fileParts=0;
                
-               //fprintf(stderr,"FSPARTICLE debug set %s , tot%d mask=%d \n", filename, totpart, readMask      );
-               for(a=0; a<totpart; a++, ftime+=dtime) {
+               for(a=0; a<totpart; a++) {
                        int ptype=0;
                        short shsize=0;
                        float convertSize=0.0;
                        gzread(gzf, &ptype, sizeof( ptype )); 
-                       //if(a<25) fprintf(stderr,"FSPARTICLE debug set %s , a%d t=%d , mask=%d  , active%d\n", filename, a, ptype, readMask, activeParts       );
                        if(ptype&readMask) {
                                activeParts++;
                                pa= new_particle(paf);
                                pa->time= ftime;
-                               pa->lifetime= ftime + G.scene->r.efra +1.0;
+                               pa->lifetime= ftime + 10000.; // add large number to make sure they are displayed, G.scene->r.efra +1.0;
                                pa->co[0] = 0.0;
                                pa->co[1] = 
                                pa->co[2] = 1.0*(float)a / (float)totpart;
@@ -1692,18 +1693,20 @@ void build_particle_system(Object *ob)
                                        gzread(gzf, &wrf, sizeof( wrf )); 
                                        vel[j] = wrf;
                                }
+                               //if(a<25) fprintf(stderr,"FSPARTICLE debug set %s , a%d = %f,%f,%f , life=%f \n", filename, a, pa->co[0],pa->co[1],pa->co[2], pa->lifetime );
                        } else {
                                // skip...
                                for(j=0; j<2*3+1; j++) {
                                        float wrf; gzread(gzf, &wrf, sizeof( wrf )); 
                                }
                        }
-                       //fprintf(stderr,"FSPARTICLE debug set %s , a%d = %f,%f,%f , life=%f \n", filename, a, pa->co[0],pa->co[1],pa->co[2], pa->lifetime );
+                       fileParts++;
                }
                gzclose( gzf );
 
                totpart = paf->totpart = activeParts;
-               //fprintf(stderr,"PARTOBH debug  %s %d \n", ob->id.name, totpart); // DEBUG
+               snprintf(debugStrBuffer,256,"readFsPartData::done - particles:%d, active:%d, file:%d, mask:%d  \n", paf->totpart,activeParts,fileParts,readMask);
+               elbeemDebugOut(debugStrBuffer);
                return;
        }
        
index 2740c7c..6427fc1 100644 (file)
@@ -2531,6 +2531,7 @@ static void direct_link_object(FileData *fd, Object *ob)
                ob->fluidsimSettings->orgMesh = NULL; //ob->data;
                ob->fluidsimSettings->meshSurface = NULL;
                ob->fluidsimSettings->meshBB = NULL;
+               ob->fluidsimSettings->meshSurfNormals = NULL;
        }
        
        link_list(fd, &ob->prop);
index e2e2ecd..eb81ff7 100644 (file)
@@ -96,14 +96,19 @@ typedef struct FluidsimSettings {
 
        /* additional flags depending on the type, lower short contains flags
         * to check validity, higher short additional flags */
-       int typeFlags;
+       short typeFlags;
+       char  domainNovecgen,dummyc;
 
        /* boundary "stickiness" for part slip values */
        float partSlipValue;
        /* particle generation - on if >0, then determines amount */
-       float generateParticles, dummy;
+       float generateParticles;
+       /* smooth fluid surface? */
+       float surfaceSmoothing;
        /* particle display - size scaling, and alpha influence */
        float particleInfSize, particleInfAlpha;
+       /* testing vars */
+       float farFieldSize,dummyf;
 
        /* save fluidsurface normals in mvert.no, and surface vertex velocities (if available) in mvert.co */
        struct MVert *meshSurfNormals;
@@ -119,13 +124,12 @@ typedef struct FluidsimSettings {
 #define OB_FLUIDSIM_OUTFLOW     32
 #define OB_FLUIDSIM_PARTICLE    64
 
-#define OB_TYPEFLAG_START       16
+#define OB_TYPEFLAG_START       0
 #define OB_FSGEO_THIN           (1<<(OB_TYPEFLAG_START+1))
 #define OB_FSBND_NOSLIP         (1<<(OB_TYPEFLAG_START+2))
 #define OB_FSBND_PARTSLIP       (1<<(OB_TYPEFLAG_START+3))
 #define OB_FSBND_FREESLIP       (1<<(OB_TYPEFLAG_START+4))
 #define OB_FSINFLOW_LOCALCOORD  (1<<(OB_TYPEFLAG_START+5))
-#define OB_FSDOMAIN_NOVECGEN      (1<<(OB_TYPEFLAG_START+6))
 
 // guiDisplayMode particle flags
 #define OB_FSDOM_GEOM     1
index afd1056..3f11598 100644 (file)
@@ -927,7 +927,7 @@ static void render_particle_system(Render *re, Object *ob, PartEff *paf)
                if(useFluidsimParticles) {
                        // rescale to 1.0-10.0, then div by 5 afterwards, gives values in range 0.2-2.0
                        double fspsize = ((double)pa->rt / 1000.0f) / 5.0 ; 
-                       haloScale = (float)pow(fspsize, (double)ob->fluidsimSettings->particleInfSize);
+                       haloScale = 1.0/(float)pow(fspsize, (double)ob->fluidsimSettings->particleInfSize);
                        ma->alpha = iniAlpha / (float)pow( fspsize, (double)ob->fluidsimSettings->particleInfAlpha);
                        if(ma->alpha>1.) ma->alpha = 1.;
                }
index 1179d7c..adba735 100644 (file)
@@ -2454,25 +2454,31 @@ static void object_panel_fluidsim(Object *ob)
                                        uiDefBut(block, LABEL, 0, "Domain boundary type settings:",             0,yline,300,objHeight, NULL, 0.0, 0, 0, 0, "");
                                        yline -= lineHeight;
                                        uiBlockBeginAlign(block);
-                                       uiDefButI(block, ROW, REDRAWBUTSOBJECT ,"Noslip",   00, yline,100,objHeight, &fss->typeFlags, 15.0, OB_FSBND_NOSLIP,   20.0, 1.0, "Obstacle causes zero normal and tangential velocity (=sticky). Default for all. Only option for moving objects.");
-                                       uiDefButI(block, ROW, REDRAWBUTSOBJECT ,"Part",    100, yline,100,objHeight, &fss->typeFlags, 15.0, OB_FSBND_PARTSLIP, 20.0, 2.0, "Mix between no-slip and free-slip. Non moving objects only!");
-                                       uiDefButI(block, ROW, REDRAWBUTSOBJECT ,"Free",          200, yline,100,objHeight, &fss->typeFlags, 15.0, OB_FSBND_FREESLIP, 20.0, 3.0, "Obstacle only causes zero normal velocity (=not sticky). Non moving objects only!");
+                                       uiDefButS(block, ROW, REDRAWBUTSOBJECT ,"Noslip",    0, yline,100,objHeight, &fss->typeFlags, 15.0, OB_FSBND_NOSLIP,   20.0, 1.0, "Obstacle causes zero normal and tangential velocity (=sticky). Default for all. Only option for moving objects.");
+                                       uiDefButS(block, ROW, REDRAWBUTSOBJECT ,"Part",    100, yline,100,objHeight, &fss->typeFlags, 15.0, OB_FSBND_PARTSLIP, 20.0, 2.0, "Mix between no-slip and free-slip. Non moving objects only!");
+                                       uiDefButS(block, ROW, REDRAWBUTSOBJECT ,"Free",          200, yline,100,objHeight, &fss->typeFlags, 15.0, OB_FSBND_FREESLIP, 20.0, 3.0, "Obstacle only causes zero normal velocity (=not sticky). Non moving objects only!");
                                        uiBlockEndAlign(block);
                                        yline -= lineHeight;
+
+                                       uiDefBut(block, LABEL, 0, "PartSlipValue:",             0,yline,150,objHeight, NULL, 0.0, 0, 0, 0, "");
                                        if(fss->typeFlags&OB_FSBND_PARTSLIP) {
-                                               uiDefBut(block, LABEL, 0, "PartSlipValue:",             0,yline,150,objHeight, NULL, 0.0, 0, 0, 0, "");
                                                uiDefButF(block, NUM, B_DIFF, "", 150, yline,150,objHeight, &fss->partSlipValue, 0.0, 0.1, 10,0, ".");
-                                               yline -= lineHeight;
-                                       }
+                                       } else { uiDefBut(block, LABEL, 0, "-", 150,yline,150,objHeight, NULL, 0.0, 0, 0, 0, ""); }
+                                       yline -= lineHeight;
+                                       // copied from obstacle...
 
                                        uiDefBut(block, LABEL, 0, "Generate Particles:",                0,yline,200,objHeight, NULL, 0.0, 0, 0, 0, "");
                                        uiDefButF(block, NUM, B_DIFF, "", 200, yline,100,objHeight, &fss->generateParticles, 0.0, 10.0, 10,0, "Amount of particles to generate (0=off, 1=normal, >1=more).");
                                        yline -= lineHeight;
 
+                                       uiDefBut(block, LABEL, 0, "Surface Smoothing:",         0,yline,200,objHeight, NULL, 0.0, 0, 0, 0, "");
+                                       uiDefButF(block, NUM, B_DIFF, "", 200, yline,100,objHeight, &fss->surfaceSmoothing, 0.0, 5.0, 10,0, "Amount of surface smoothing (0=off, 1=normal, >1=stronger smoothing).");
+                                       yline -= lineHeight;
+
+                                       // use new variable...
                                        uiDefBut(block, LABEL, 0, "Generate&Use SpeedVecs:",            0,yline,200,objHeight, NULL, 0.0, 0, 0, 0, "");
-                                 uiDefButBitI(block, TOG, OB_FSDOMAIN_NOVECGEN, REDRAWBUTSOBJECT, "Disable",     200, yline,100,objHeight, &fss->typeFlags, 0, 0, 0, 0, "Default is to generate and use fluidsim vertex speed vectors, this option switches calculation off during bake, and disables loading.");
+                                 uiDefButBitC(block, TOG, 1, REDRAWBUTSOBJECT, "Disable",     200, yline,100,objHeight, &fss->domainNovecgen, 0, 0, 0, 0, "Default is to generate and use fluidsim vertex speed vectors, this option switches calculation off during bake, and disables loading.");
                                        yline -= lineHeight;
-                                       // copied from obstacle...
                                }
                        }
                        else if(
@@ -2492,7 +2498,7 @@ static void object_panel_fluidsim(Object *ob)
 
                                if(fss->type == OB_FLUIDSIM_INFLOW) {
                                        uiDefBut(block, LABEL, 0, "Local Inflow Coords",                0,yline,200,objHeight, NULL, 0.0, 0, 0, 0, "");
-                                 uiDefButBitI(block, TOG, OB_FSINFLOW_LOCALCOORD, REDRAWBUTSOBJECT, "Enable",     200, yline,100,objHeight, &fss->typeFlags, 0, 0, 0, 0, "Use local coordinates for inflow (e.g. for rotating objects).");
+                                 uiDefButBitS(block, TOG, OB_FSINFLOW_LOCALCOORD, REDRAWBUTSOBJECT, "Enable",     200, yline,100,objHeight, &fss->typeFlags, 0, 0, 0, 0, "Use local coordinates for inflow (e.g. for rotating objects).");
                                  yline -= lineHeight;
                                }
                        }
@@ -2504,31 +2510,26 @@ static void object_panel_fluidsim(Object *ob)
                                yline -= lineHeight + 5;
 
                                uiBlockBeginAlign(block);
-                               uiDefButI(block, ROW, REDRAWBUTSOBJECT ,"Noslip",   00, yline,100,objHeight, &fss->typeFlags, 15.0, OB_FSBND_NOSLIP,   20.0, 1.0, "Obstacle causes zero normal and tangential velocity (=sticky). Default for all. Only option for moving objects.");
-                               uiDefButI(block, ROW, REDRAWBUTSOBJECT ,"Part",    100, yline,100,objHeight, &fss->typeFlags, 15.0, OB_FSBND_PARTSLIP, 20.0, 2.0, "Mix between no-slip and free-slip. Non moving objects only!");
-                               uiDefButI(block, ROW, REDRAWBUTSOBJECT ,"Free",          200, yline,100,objHeight, &fss->typeFlags, 15.0, OB_FSBND_FREESLIP, 20.0, 3.0, "Obstacle only causes zero normal velocity (=not sticky). Non moving objects only!");
+                               uiDefButS(block, ROW, REDRAWBUTSOBJECT ,"Noslip",   00, yline,100,objHeight, &fss->typeFlags, 15.0, OB_FSBND_NOSLIP,   20.0, 1.0, "Obstacle causes zero normal and tangential velocity (=sticky). Default for all. Only option for moving objects.");
+                               uiDefButS(block, ROW, REDRAWBUTSOBJECT ,"Part",    100, yline,100,objHeight, &fss->typeFlags, 15.0, OB_FSBND_PARTSLIP, 20.0, 2.0, "Mix between no-slip and free-slip. Non moving objects only!");
+                               uiDefButS(block, ROW, REDRAWBUTSOBJECT ,"Free",          200, yline,100,objHeight, &fss->typeFlags, 15.0, OB_FSBND_FREESLIP, 20.0, 3.0, "Obstacle only causes zero normal velocity (=not sticky). Non moving objects only!");
                                uiBlockEndAlign(block);
                                yline -= lineHeight;
 
+                               uiDefBut(block, LABEL, 0, "PartSlipValue:",             0,yline,150,objHeight, NULL, 0.0, 0, 0, 0, "");
                                if(fss->typeFlags&OB_FSBND_PARTSLIP) {
-                                       uiDefBut(block, LABEL, 0, "PartSlipValue:",             0,yline,150,objHeight, NULL, 0.0, 0, 0, 0, "");
                                        uiDefButF(block, NUM, B_DIFF, "", 150, yline,150,objHeight, &fss->partSlipValue, 0.0, 0.1, 10,0, ".");
-                                       yline -= lineHeight;
-                               }
+                               } else { uiDefBut(block, LABEL, 0, "-", 150,yline,150,objHeight, NULL, 0.0, 0, 0, 0, ""); }
+                               yline -= lineHeight;
 
                                yline -= lineHeight;
                        }
                        else if(fss->type == OB_FLUIDSIM_PARTICLE) {
 
-                               if(fss->guiDisplayMode==0) fss->guiDisplayMode=2; // default drops
-                               uiDefBut(block, LABEL,   0, "Part.-Type:",               0,yline,100,objHeight, NULL, 0.0, 0, 0, 0, "");
-                               // TODO make toggle buttons
-                               //uiDefButS(block, MENU, B_FLUIDSIM_FORCEREDRAW, "Gui%t|Bubble %x2|Drop %x4|Newparts %x8|Float %x16",   
-                                                //100,yline,200,objHeight, &fss->guiDisplayMode, 0, 0, 0, 0, "Which type of particles to display.");
-                               //uiDefButS(block, MENU, B_DIFF, "Render%t|Geometry %x1|Preview %x2|Final %x3", 
-                                               //195,yline,105,objHeight, &fss->renderDisplayMode, 0, 0, 0, 0, "How to display the fluid mesh for rendering.");
-                               uiDefBut(block, LABEL,   0, "Drops",    100,yline,200,objHeight, NULL, 0.0, 0, 0, 0, "");
-                               fss->guiDisplayMode = 4; // fix to drops for now
+                               if(fss->typeFlags==0) fss->typeFlags=4; // default drops
+                               fss->typeFlags = 4; // fix to drops for now
+                               uiDefBut(block, LABEL,   0, "Part.-Type:",    0,yline,100,objHeight, NULL, 0.0, 0, 0, 0, "");
+                               uiDefBut(block, LABEL,   0, "Drops",        100,yline,200,objHeight, NULL, 0.0, 0, 0, 0, "");
                                yline -= lineHeight;
 
                                uiDefBut(block, LABEL, 0, "Size Influence:",            0,yline,150,objHeight, NULL, 0.0, 0, 0, 0, "");
index ef1d5fc..15e41be 100644 (file)
@@ -84,6 +84,7 @@
 #include "BSE_headerbuttons.h"
 
 #include "mydevice.h"
+#include "blendef.h"
 #include "SDL.h"
 #include "SDL_thread.h"
 #include "SDL_mutex.h"
@@ -109,9 +110,6 @@ void initElbeemMesh(struct Object *ob, int *numVertices, float **vertices, int *
 extern int start_progress_bar(void);
 extern void end_progress_bar(void);
 extern int progress_bar(float done, char *busy_info);
-// global solver state
-extern int gElbeemState;
-extern char gElbeemErrorString[];
 
 double fluidsimViscosityPreset[6] = {
        -1.0,   /* unused */
@@ -197,9 +195,11 @@ FluidsimSettings *fluidsimSettingsNew(struct Object *srcob)
        fluidsimGetAxisAlignedBB(srcob->data, srcob->obmat, fss->bbStart, fss->bbSize, &fss->meshBB);
        
        fss->typeFlags = 0;
+       fss->domainNovecgen = 0;
        fss->partSlipValue = 0.0;
 
        fss->generateParticles = 0.0;
+       fss->surfaceSmoothing = 1.0;
        fss->particleInfSize = 0.0;
        fss->particleInfAlpha = 0.0;
 
@@ -445,6 +445,24 @@ void fluidsimBake(struct Object *ob)
                return;
        }
 
+       /* no object pointer, find in selected ones.. */
+       if(!ob) {
+               Base *base;
+               for(base=G.scene->base.first; base; base= base->next) {
+                       if ( ((base)->flag & SELECT) 
+                                       // ignore layer setting for now? && ((base)->lay & G.vd->lay) 
+                                ) {
+                               if((!ob)&&(base->object->fluidsimFlag & OB_FLUIDSIM_ENABLE)&&(base->object->type==OB_MESH)) {
+                                       if(base->object->fluidsimSettings->type == OB_FLUIDSIM_DOMAIN) {
+                                               ob = base->object;
+                                       }
+                               }
+                       }
+               }
+               // no domains found?
+               if(!ob) return;
+       }
+
        /* check if there's another domain... */
        for(obit= G.main->object.first; obit; obit= obit->id.next) {
                if((obit->fluidsimFlag & OB_FLUIDSIM_ENABLE)&&(obit->type==OB_MESH)) {
@@ -764,6 +782,8 @@ void fluidsimBake(struct Object *ob)
                fsset.gstar = domainSettings->gstar;
                fsset.maxRefine = domainSettings->maxRefine; // check <-> gridlevels
                fsset.generateParticles = domainSettings->generateParticles; 
+               fsset.surfaceSmoothing = domainSettings->surfaceSmoothing; 
+               fsset.farFieldSize = domainSettings->farFieldSize; 
                strcpy( fsset.outputPath, targetFile);
 
                // domain channels
@@ -778,8 +798,7 @@ void fluidsimBake(struct Object *ob)
                else if((domainSettings->typeFlags&OB_FSBND_PARTSLIP)) fsset.obstacleType = FLUIDSIM_OBSTACLE_PARTSLIP;
                else if((domainSettings->typeFlags&OB_FSBND_FREESLIP)) fsset.obstacleType = FLUIDSIM_OBSTACLE_FREESLIP;
                fsset.obstaclePartslip = domainSettings->partSlipValue;
-               fsset.generateVertexVectors = (int)(!(domainSettings->typeFlags&OB_FSDOMAIN_NOVECGEN));
-               // fprintf(stderr," VVV %d %d \n",fsset.generateVertexVectors , (domainSettings->typeFlags&OB_FSDOMAIN_NOVECGEN)); // DEBUG
+               fsset.generateVertexVectors = (domainSettings->domainNovecgen==0);
 
                // init blender trafo matrix
                // fprintf(stderr,"elbeemInit - mpTrafo:\n");
@@ -950,8 +969,8 @@ void fluidsimBake(struct Object *ob)
                        "  size = " "%d"            /* gridSize*/ ";  \n" 
                        "  surfacepreview = " "%d"  /* previewSize*/ "; \n" 
                        "  dump_velocities = " "%d"  /* vector dump */ "; \n" 
-                       "  smoothsurface = 1.0;  \n"
-                       "  smoothnormals = 1.0;  \n"
+                       "  smoothsurface = %f;  \n"  /* smoothing */
+                       "  smoothnormals = %f;  \n"
                        "  geoinitid = 1;  \n"  "\n" 
                        "  isovalue =  0.4900; \n" 
                        "  isoweightmethod = 1; \n"  "\n" ;
@@ -959,7 +978,7 @@ void fluidsimBake(struct Object *ob)
                        fprintf(fileCfg, simString,
                                        (double)domainSettings->realsize, (double)domainSettings->animStart, (double)domainSettings->gstar,
                                        gridlevels, (int)domainSettings->resolutionxyz, (int)domainSettings->previewresxyz ,
-                                       (int)(!(domainSettings->typeFlags&OB_FSDOMAIN_NOVECGEN))
+                     (int)(domainSettings->domainNovecgen==0), domainSettings->surfaceSmoothing
                                        );
 
                        if((domainSettings->typeFlags&OB_FSBND_NOSLIP)) bi=0;
@@ -972,6 +991,7 @@ void fluidsimBake(struct Object *ob)
                        fluidsimPrintChannel(fileCfg, channelDomainViscosity,allchannelSize,"p_viscosity",CHANNEL_FLOAT);
                        fluidsimPrintChannel(fileCfg, channelDomainGravity,  allchannelSize,"p_gravity",CHANNEL_VEC);
 
+                       fprintf(fileCfg, "  partgenprob = %f; \n", domainSettings->generateParticles); // debug test
                        fprintf(fileCfg,  "\n} \n" );
                }
 
@@ -1199,17 +1219,16 @@ void fluidsimBake(struct Object *ob)
 
        if(!simAborted) {
                char fsmessage[512];
+               char elbeemerr[256];
                strcpy(fsmessage,"Fluidsim Bake Error: ");
                // check if some error occurred
                if(globalBakeState==-2) {
                        strcat(fsmessage,"Failed to initialize [Msg: ");
-#ifndef WIN32
-                       // msvc seems to have problem accessing the gElbeemErrorString var
-                       strcat(fsmessage,"[Msg: ");
-                       strcat(fsmessage,gElbeemErrorString);
-                       strcat(fsmessage,"]");
-#endif // WIN32
-                       strcat(fsmessage,"|OK%x0");
+
+                       elbeemGetErrorString(elbeemerr);
+                       strcat(fsmessage,elbeemerr);
+
+                       strcat(fsmessage,"] |OK%x0");
                        pupmenu(fsmessage);
                } // init error
        }
index b16b5c7..d466984 100644 (file)
@@ -1216,7 +1216,10 @@ static void winqreadview3dspace(ScrArea *sa, void *spacedata, BWinEvent *evt)
                                else if(G.qual==LR_CTRLKEY) {
                                        if(okee("Bake all selected")) {
                                                extern void softbody_bake(Object *ob);
+                                               extern void fluidsimBake(Object *ob);
                                                softbody_bake(NULL);
+                                               // also bake first domain of selected objects...
+                                               fluidsimBake(NULL);
                                        }
                                }
                                else if(G.qual==0)