blender-archive/source/blender/editors/physics/physics_fluid.c

/*
 * ***** BEGIN GPL LICENSE BLOCK *****
 *
 * This program is free software; you can redistribute it and/or
 * modify it under the terms of the GNU General Public License
 * as published by the Free Software Foundation; either version 2
 * of the License, or (at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software Foundation,
 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
 *
 * The Original Code is Copyright (C) Blender Foundation
 * All rights reserved.
 *
 * The Original Code is: all of this file.
 *
 * Contributor(s): none yet.
 *
 * ***** END GPL LICENSE BLOCK *****
 */

/** \file blender/editors/physics/physics_fluid.c
 *  \ingroup edphys
 */

#include <math.h>
#include <stdlib.h>
#include <string.h>
#include <sys/stat.h>

#include "MEM_guardedalloc.h"

/* types */
#include "DNA_action_types.h"
#include "DNA_object_types.h"
#include "DNA_object_fluidsim_types.h"	

#include "BLI_blenlib.h"
#include "BLI_math.h"
#include "BLI_utildefines.h"

#include "BKE_context.h"
#include "BKE_customdata.h"
#include "BKE_fluidsim.h"
#include "BKE_modifier.h"
#include "BKE_object.h"
#include "BKE_report.h"
#include "BKE_scene.h"

#include "DEG_depsgraph.h"

#include "LBM_fluidsim.h"

#include "ED_screen.h"

#include "WM_types.h"
#include "WM_api.h"

#include "physics_intern.h" // own include

/* enable/disable overall compilation */
#ifdef WITH_MOD_FLUID

#include "BKE_global.h"

#include "WM_api.h"

#include "DNA_scene_types.h"
#include "DNA_mesh_types.h"


static float get_fluid_viscosity(FluidsimSettings *settings)
{
	return (1.0f/powf(10.0f, settings->viscosityExponent)) * settings->viscosityValue;
}

static float get_fluid_rate(FluidsimSettings *settings)
{
	float rate = 1.0f; /* default rate if not animated... */
	
	rate = settings->animRate;
	
	if (rate < 0.0f)
		rate = 0.0f;
	
	return rate;
}

static void get_fluid_gravity(float *gravity, Scene *scene, FluidsimSettings *fss)
{
	if (scene->physics_settings.flag & PHYS_GLOBAL_GRAVITY) {
		copy_v3_v3(gravity, scene->physics_settings.gravity);
	}
	else {
		copy_v3_v3(gravity, fss->grav);
	}
}

static float get_fluid_size_m(Scene *scene, Object *domainob, FluidsimSettings *fss)
{
	if (!scene->unit.system) {
		return fss->realsize;
	}
	else {
		float dim[3];
		float longest_axis;
		
		BKE_object_dimensions_get(domainob, dim);
		longest_axis = max_fff(dim[0], dim[1], dim[2]);
		
		return longest_axis * scene->unit.scale_length;
	}
}

static bool fluid_is_animated_mesh(FluidsimSettings *fss)
{
	return ((fss->type == OB_FLUIDSIM_CONTROL) || fss->domainNovecgen);
}

/* ********************** fluid sim settings struct functions ********************** */

#if 0
/* helper function */
void fluidsimGetGeometryObjFilename(Object *ob, char *dst)  //, char *srcname)
{
	//BLI_snprintf(dst, FILE_MAXFILE, "%s_cfgdata_%s.bobj.gz", srcname, ob->id.name);
	BLI_snprintf(dst, FILE_MAXFILE, "fluidcfgdata_%s.bobj.gz", ob->id.name);
}
#endif


/* ********************** fluid sim channel helper functions ********************** */

typedef struct FluidAnimChannels {
	int length;
	
	double aniFrameTime;
	
	float *timeAtFrame;
	float *DomainTime;
	float *DomainGravity;
	float *DomainViscosity;
} FluidAnimChannels;

typedef struct FluidObject {
	struct FluidObject *next, *prev;
	
	struct Object *object;
	
	float *Translation;
	float *Rotation;
	float *Scale;
	float *Active;
	
	float *InitialVelocity;
	
	float *AttractforceStrength;
	float *AttractforceRadius;
	float *VelocityforceStrength;
	float *VelocityforceRadius;
	
	float *VertexCache;
	int numVerts, numTris;
} FluidObject;

// no. of entries for the two channel sizes
#define CHANNEL_FLOAT 1
#define CHANNEL_VEC   3

// simplify channels before printing
// for API this is done anyway upon init
#if 0
static void fluidsimPrintChannel(FILE *file, float *channel, int paramsize, char *str, int entries) 
{ 
	int i, j;
	int channelSize = paramsize; 

	if (entries == 3) {
		elbeemSimplifyChannelVec3(channel, &channelSize);
	}
	else if (entries == 1) {
		elbeemSimplifyChannelFloat(channel, &channelSize);
	}
	else {
		/* invalid, cant happen? */
	}

	fprintf(file, "      CHANNEL %s =\n", str);
	for (i=0; i < channelSize; i++) {
		fprintf(file, "        ");
		for (j=0;j <= entries;j++) {  // also print time value
			fprintf(file, " %f ", channel[i*(entries + 1) + j]);
			if (j == entries-1) { fprintf(file, "  "); }
		}
		fprintf(file, "\n");
	}

	fprintf(file,  "      ;\n");
}
#endif


/* Note: fluid anim channel data layout
 * ------------------------------------
 * CHANNEL_FLOAT:
 * frame 1     |frame 2
 * [dataF][time][dataF][time]
 *
 * CHANNEL_VEC:
 * frame 1                   |frame 2
 * [dataX][dataY][dataZ][time][dataX][dataY][dataZ][time]
 *
 */

static void init_time(FluidsimSettings *domainSettings, FluidAnimChannels *channels)
{
	int i;
	
	channels->timeAtFrame = MEM_callocN((channels->length + 1) * sizeof(float), "timeAtFrame channel");
	
	channels->timeAtFrame[0] = channels->timeAtFrame[1] = domainSettings->animStart; // start at index 1
	
	for (i=2; i <= channels->length; i++) {
		channels->timeAtFrame[i] = channels->timeAtFrame[i - 1] + (float)channels->aniFrameTime;
	}
}

/* if this is slow, can replace with faster, less readable code */
static void set_channel(float *channel, float time, float *value, int i, int size)
{
	if (size == CHANNEL_FLOAT) {
		channel[(i * 2) + 0] = value[0];
		channel[(i * 2) + 1] = time;
	}
	else if (size == CHANNEL_VEC) {
		channel[(i * 4) + 0] = value[0];
		channel[(i * 4) + 1] = value[1];
		channel[(i * 4) + 2] = value[2];
		channel[(i * 4) + 3] = time;
	}
}

static void set_vertex_channel(Depsgraph *depsgraph, float *channel, float time, struct Scene *scene, struct FluidObject *fobj, int i)
{
	Object *ob = fobj->object;
	FluidsimModifierData *fluidmd = (FluidsimModifierData *)modifiers_findByType(ob, eModifierType_Fluidsim);
	float *verts;
	int *tris=NULL, numVerts=0, numTris=0;
	int modifierIndex = BLI_findindex(&ob->modifiers, fluidmd);
	int framesize = (3*fobj->numVerts) + 1;
	int j;
	
	if (channel == NULL)
		return;
	
	initElbeemMesh(depsgraph, scene, ob, &numVerts, &verts, &numTris, &tris, 1, modifierIndex);
	
	/* don't allow mesh to change number of verts in anim sequence */
	if (numVerts != fobj->numVerts) {
		MEM_freeN(channel);
		channel = NULL;
		return;
	}
	
	/* fill frame of channel with vertex locations */
	for (j=0; j < (3*numVerts); j++) {
		channel[i*framesize + j] = verts[j];
	}
	channel[i*framesize + framesize-1] = time;
	
	MEM_freeN(verts);
	MEM_freeN(tris);
}

static void free_domain_channels(FluidAnimChannels *channels)
{
	if (!channels->timeAtFrame)
		return;
	MEM_freeN(channels->timeAtFrame);
	channels->timeAtFrame = NULL;
	MEM_freeN(channels->DomainGravity);
	channels->DomainGravity = NULL;
	MEM_freeN(channels->DomainViscosity);
	channels->DomainViscosity = NULL;
	MEM_freeN(channels->DomainTime);
	channels->DomainTime = NULL;
}

static void free_all_fluidobject_channels(ListBase *fobjects)
{
	FluidObject *fobj;
	
	for (fobj=fobjects->first; fobj; fobj=fobj->next) {
		if (fobj->Translation) {
			MEM_freeN(fobj->Translation);
			fobj->Translation = NULL;
			MEM_freeN(fobj->Rotation);
			fobj->Rotation = NULL;
			MEM_freeN(fobj->Scale);
			fobj->Scale = NULL;
			MEM_freeN(fobj->Active);
			fobj->Active = NULL;
			MEM_freeN(fobj->InitialVelocity);
			fobj->InitialVelocity = NULL;
		}
		
		if (fobj->AttractforceStrength) {
			MEM_freeN(fobj->AttractforceStrength);
			fobj->AttractforceStrength = NULL;
			MEM_freeN(fobj->AttractforceRadius);
			fobj->AttractforceRadius = NULL;
			MEM_freeN(fobj->VelocityforceStrength);
			fobj->VelocityforceStrength = NULL;
			MEM_freeN(fobj->VelocityforceRadius);
			fobj->VelocityforceRadius = NULL;
		}
		
		if (fobj->VertexCache) {
			MEM_freeN(fobj->VertexCache);
			fobj->VertexCache = NULL;
		}
	}
}

static void fluid_init_all_channels(bContext *C, Object *UNUSED(fsDomain), FluidsimSettings *domainSettings, FluidAnimChannels *channels, ListBase *fobjects)
{
	Scene *scene = CTX_data_scene(C);
	ViewLayer *view_layer = CTX_data_view_layer(C);
	Depsgraph *depsgraph = CTX_data_depsgraph(C);
	Base *base;
	int i;
	int length = channels->length;
	float eval_time;

	/* init time values (assuming that time moves at a constant speed; may be overridden later) */
	init_time(domainSettings, channels);
	
	/* allocate domain animation channels */
	channels->DomainGravity = MEM_callocN(length * (CHANNEL_VEC+1) * sizeof(float), "channel DomainGravity");
	channels->DomainViscosity = MEM_callocN(length * (CHANNEL_FLOAT+1) * sizeof(float), "channel DomainViscosity");
	channels->DomainTime = MEM_callocN(length * (CHANNEL_FLOAT+1) * sizeof(float), "channel DomainTime");
	
	/* allocate fluid objects */
	for (base = FIRSTBASE(view_layer); base; base = base->next) {
		Object *ob = base->object;
		FluidsimModifierData *fluidmd = (FluidsimModifierData *)modifiers_findByType(ob, eModifierType_Fluidsim);
		
		if (fluidmd) {
			FluidObject *fobj = MEM_callocN(sizeof(FluidObject), "Fluid Object");
			fobj->object = ob;
			
			if (ELEM(fluidmd->fss->type, OB_FLUIDSIM_DOMAIN, OB_FLUIDSIM_PARTICLE)) {
				BLI_addtail(fobjects, fobj);
				continue;
			}
			
			fobj->Translation = MEM_callocN(length * (CHANNEL_VEC+1) * sizeof(float), "fluidobject Translation");
			fobj->Rotation = MEM_callocN(length * (CHANNEL_VEC+1) * sizeof(float), "fluidobject Rotation");
			fobj->Scale = MEM_callocN(length * (CHANNEL_VEC+1) * sizeof(float), "fluidobject Scale");
			fobj->Active = MEM_callocN(length * (CHANNEL_FLOAT+1) * sizeof(float), "fluidobject Active");
			fobj->InitialVelocity = MEM_callocN(length * (CHANNEL_VEC+1) * sizeof(float), "fluidobject InitialVelocity");
			
			if (fluidmd->fss->type == OB_FLUIDSIM_CONTROL) {
				fobj->AttractforceStrength = MEM_callocN(length * (CHANNEL_FLOAT+1) * sizeof(float), "fluidobject AttractforceStrength");
				fobj->AttractforceRadius = MEM_callocN(length * (CHANNEL_FLOAT+1) * sizeof(float), "fluidobject AttractforceRadius");
				fobj->VelocityforceStrength = MEM_callocN(length * (CHANNEL_FLOAT+1) * sizeof(float), "fluidobject VelocityforceStrength");
				fobj->VelocityforceRadius = MEM_callocN(length * (CHANNEL_FLOAT+1) * sizeof(float), "fluidobject VelocityforceRadius");
			}
			
			if (fluid_is_animated_mesh(fluidmd->fss)) {
				float *verts=NULL;
				int *tris=NULL, modifierIndex = BLI_findindex(&ob->modifiers, (ModifierData *)fluidmd);

				initElbeemMesh(depsgraph, scene, ob, &fobj->numVerts, &verts, &fobj->numTris, &tris, 0, modifierIndex);
				fobj->VertexCache = MEM_callocN(length *((fobj->numVerts*CHANNEL_VEC)+1) * sizeof(float), "fluidobject VertexCache");
				
				MEM_freeN(verts);
				MEM_freeN(tris);
			}
			
			BLI_addtail(fobjects, fobj);
		}
	}
	
	/* now we loop over the frames and fill the allocated channels with data */
	for (i=0; i < channels->length; i++) {
		FluidObject *fobj;
		float viscosity, gravity[3];
		float timeAtFrame, time;
		
		eval_time = domainSettings->bakeStart + i;
		
		/* XXX: This can't be used due to an anim sys optimization that ignores recalc object animation,
		 * leaving it for the depgraph (this ignores object animation such as modifier properties though... :/ )
		 * --> BKE_animsys_evaluate_all_animation(G.main, eval_time);
		 * This doesn't work with drivers:
		 * --> BKE_animsys_evaluate_animdata(&fsDomain->id, fsDomain->adt, eval_time, ADT_RECALC_ALL);
		 */
		
		/* Modifying the global scene isn't nice, but we can do it in 
		 * this part of the process before a threaded job is created */
		scene->r.cfra = (int)eval_time;
		ED_update_for_newframe(CTX_data_main(C), depsgraph);
		
		/* now scene data should be current according to animation system, so we fill the channels */
		
		/* Domain time */
		// TODO: have option for not running sim, time mangling, in which case second case comes in handy
		if (channels->DomainTime) {
			time = get_fluid_rate(domainSettings) * (float)channels->aniFrameTime;
			timeAtFrame = channels->timeAtFrame[i] + time;
			
			channels->timeAtFrame[i+1] = timeAtFrame;
			set_channel(channels->DomainTime, i, &time, i, CHANNEL_FLOAT);
		}
		else {
			timeAtFrame = channels->timeAtFrame[i+1];
		}
		
		/* Domain properties - gravity/viscosity */
		get_fluid_gravity(gravity, scene, domainSettings);
		set_channel(channels->DomainGravity, timeAtFrame, gravity, i, CHANNEL_VEC);
		viscosity = get_fluid_viscosity(domainSettings);
		set_channel(channels->DomainViscosity, timeAtFrame, &viscosity, i, CHANNEL_FLOAT);
		
		/* object movement */
		for (fobj=fobjects->first; fobj; fobj=fobj->next) {
			Object *ob = fobj->object;
			FluidsimModifierData *fluidmd = (FluidsimModifierData *)modifiers_findByType(ob, eModifierType_Fluidsim);
			float active= (float) ((fluidmd->fss->flag & OB_FLUIDSIM_ACTIVE) > 0 ? 1 : 0);
			float rot_d[3] = {0.f, 0.f, 0.f}, old_rot[3] = {0.f, 0.f, 0.f};
			
			if (ELEM(fluidmd->fss->type, OB_FLUIDSIM_DOMAIN, OB_FLUIDSIM_PARTICLE))
				continue;
			
			/* init euler rotation values and convert to elbeem format */
			/* get the rotation from ob->obmat rather than ob->rot to account for parent animations */
			if (i) {
				copy_v3_v3(old_rot, fobj->Rotation + 4*(i-1));
				mul_v3_fl(old_rot, (float)-M_PI / 180.f);
			}

			mat4_to_compatible_eulO(rot_d, old_rot, 0, ob->obmat);
			mul_v3_fl(rot_d, -180.0f / (float)M_PI);
			
			set_channel(fobj->Translation, timeAtFrame, ob->loc, i, CHANNEL_VEC);
			set_channel(fobj->Rotation, timeAtFrame, rot_d, i, CHANNEL_VEC);
			set_channel(fobj->Scale, timeAtFrame, ob->size, i, CHANNEL_VEC);
			set_channel(fobj->Active, timeAtFrame, &active, i, CHANNEL_FLOAT);
			set_channel(fobj->InitialVelocity, timeAtFrame, &fluidmd->fss->iniVelx, i, CHANNEL_VEC);

			// printf("Active: %f, Frame: %f\n", active, timeAtFrame);
			
			if (fluidmd->fss->type == OB_FLUIDSIM_CONTROL) {
				set_channel(fobj->AttractforceStrength, timeAtFrame, &fluidmd->fss->attractforceStrength, i, CHANNEL_FLOAT);
				set_channel(fobj->AttractforceRadius, timeAtFrame, &fluidmd->fss->attractforceRadius, i, CHANNEL_FLOAT);
				set_channel(fobj->VelocityforceStrength, timeAtFrame, &fluidmd->fss->velocityforceStrength, i, CHANNEL_FLOAT);
				set_channel(fobj->VelocityforceRadius, timeAtFrame, &fluidmd->fss->velocityforceRadius, i, CHANNEL_FLOAT);
			}
			
			if (fluid_is_animated_mesh(fluidmd->fss)) {
				set_vertex_channel(depsgraph, fobj->VertexCache, timeAtFrame, scene, fobj, i);
			}
		}
	}
}

static void export_fluid_objects(const bContext *C, ListBase *fobjects, Scene *scene, int length)
{
	Depsgraph *depsgraph = CTX_data_depsgraph(C);
	FluidObject *fobj;
	
	for (fobj=fobjects->first; fobj; fobj=fobj->next) {
		Object *ob = fobj->object;
		FluidsimModifierData *fluidmd = (FluidsimModifierData *)modifiers_findByType(ob, eModifierType_Fluidsim);
		int modifierIndex = BLI_findindex(&ob->modifiers, fluidmd);
		
		float *verts=NULL;
		int *tris=NULL;
		int numVerts=0, numTris=0;
		bool deform = fluid_is_animated_mesh(fluidmd->fss);
		
		elbeemMesh fsmesh;
		
		if (ELEM(fluidmd->fss->type, OB_FLUIDSIM_DOMAIN, OB_FLUIDSIM_PARTICLE))
			continue;
		
		elbeemResetMesh(&fsmesh);
		
		fsmesh.type = fluidmd->fss->type;
		fsmesh.name = ob->id.name;
		
		initElbeemMesh(depsgraph, scene, ob, &numVerts, &verts, &numTris, &tris, 0, modifierIndex);
		
		fsmesh.numVertices   = numVerts;
		fsmesh.numTriangles  = numTris;
		fsmesh.vertices      = verts;
		fsmesh.triangles     = tris;
		
		fsmesh.channelSizeTranslation  = 
		fsmesh.channelSizeRotation     = 
		fsmesh.channelSizeScale        = 
		fsmesh.channelSizeInitialVel   = 
		fsmesh.channelSizeActive       = length;
		
		fsmesh.channelTranslation      = fobj->Translation;
		fsmesh.channelRotation         = fobj->Rotation;
		fsmesh.channelScale            = fobj->Scale;
		fsmesh.channelActive           = fobj->Active;
		
		if ( ELEM(fsmesh.type, OB_FLUIDSIM_FLUID, OB_FLUIDSIM_INFLOW)) {
			fsmesh.channelInitialVel = fobj->InitialVelocity;
			fsmesh.localInivelCoords = ((fluidmd->fss->typeFlags & OB_FSINFLOW_LOCALCOORD) ? 1 : 0);
		}
		
		if (fluidmd->fss->typeFlags & OB_FSBND_NOSLIP)
			fsmesh.obstacleType = FLUIDSIM_OBSTACLE_NOSLIP;
		else if (fluidmd->fss->typeFlags & OB_FSBND_PARTSLIP)
			fsmesh.obstacleType = FLUIDSIM_OBSTACLE_PARTSLIP;
		else if (fluidmd->fss->typeFlags & OB_FSBND_FREESLIP)
			fsmesh.obstacleType = FLUIDSIM_OBSTACLE_FREESLIP;
		
		fsmesh.obstaclePartslip = fluidmd->fss->partSlipValue;
		fsmesh.volumeInitType = fluidmd->fss->volumeInitType;
		fsmesh.obstacleImpactFactor = fluidmd->fss->surfaceSmoothing; // misused value
		
		if (fsmesh.type == OB_FLUIDSIM_CONTROL) {
			fsmesh.cpsTimeStart = fluidmd->fss->cpsTimeStart;
			fsmesh.cpsTimeEnd = fluidmd->fss->cpsTimeEnd;
			fsmesh.cpsQuality = fluidmd->fss->cpsQuality;
			fsmesh.obstacleType = (fluidmd->fss->flag & OB_FLUIDSIM_REVERSE);
			
			fsmesh.channelSizeAttractforceRadius = 
			fsmesh.channelSizeVelocityforceStrength = 
			fsmesh.channelSizeVelocityforceRadius = 
			fsmesh.channelSizeAttractforceStrength = length;
			
			fsmesh.channelAttractforceStrength = fobj->AttractforceStrength;
			fsmesh.channelAttractforceRadius = fobj->AttractforceRadius;
			fsmesh.channelVelocityforceStrength = fobj->VelocityforceStrength;
			fsmesh.channelVelocityforceRadius = fobj->VelocityforceRadius;
		}
		else {
			fsmesh.channelAttractforceStrength =
			fsmesh.channelAttractforceRadius = 
			fsmesh.channelVelocityforceStrength = 
			fsmesh.channelVelocityforceRadius = NULL; 
		}
		
		/* animated meshes */
		if (deform) {
			fsmesh.channelSizeVertices = length;
			fsmesh.channelVertices = fobj->VertexCache;
			
			/* remove channels */
			fsmesh.channelTranslation      = 
			fsmesh.channelRotation         = 
			fsmesh.channelScale            = NULL;
			
			/* Override user settings, only noslip is supported here! */
			if (fsmesh.type != OB_FLUIDSIM_CONTROL)
				fsmesh.obstacleType = FLUIDSIM_OBSTACLE_NOSLIP;
		}
		
		elbeemAddMesh(&fsmesh);
		
		if (verts) MEM_freeN(verts);
		if (tris) MEM_freeN(tris);
	}
}

static int fluid_validate_scene(ReportList *reports, ViewLayer *view_layer, Object *fsDomain)
{
	Base *base;
	Object *newdomain = NULL;
	int channelObjCount = 0;
	int fluidInputCount = 0;

	for (base = FIRSTBASE(view_layer); base; base = base->next) {
		Object *ob = base->object;
		FluidsimModifierData *fluidmdtmp = (FluidsimModifierData *)modifiers_findByType(ob, eModifierType_Fluidsim);

		/* only find objects with fluid modifiers */
		if (!fluidmdtmp || ob->type != OB_MESH) continue;
			
		if (fluidmdtmp->fss->type == OB_FLUIDSIM_DOMAIN) {
			/* if no initial domain object given, find another potential domain */
			if (!fsDomain) {
				newdomain = ob;
			}
			/* if there's more than one domain, cancel */
			else if (fsDomain && ob != fsDomain) {
				BKE_report(reports, RPT_ERROR, "There should be only one domain object");
				return 0;
			}
		}
		
		/* count number of objects needed for animation channels */
		if ( !ELEM(fluidmdtmp->fss->type, OB_FLUIDSIM_DOMAIN, OB_FLUIDSIM_PARTICLE) )
			channelObjCount++;
		
		/* count number of fluid input objects */
		if (ELEM(fluidmdtmp->fss->type, OB_FLUIDSIM_FLUID, OB_FLUIDSIM_INFLOW))
			fluidInputCount++;
	}

	if (newdomain)
		fsDomain = newdomain;
	
	if (!fsDomain) {
		BKE_report(reports, RPT_ERROR, "No domain object found");
		return 0;
	}
	
	if (channelObjCount >= 255) {
		BKE_report(reports, RPT_ERROR, "Cannot bake with more than 256 objects");
		return 0;
	}
	
	if (fluidInputCount == 0) {
		BKE_report(reports, RPT_ERROR, "No fluid input objects in the scene");
		return 0;
	}
	
	return 1;
}


#define FLUID_SUFFIX_CONFIG		"fluidsim.cfg"
#define FLUID_SUFFIX_CONFIG_TMP	(FLUID_SUFFIX_CONFIG ".tmp")
#define FLUID_SUFFIX_SURFACE	"fluidsurface"

static bool fluid_init_filepaths(
        ReportList *reports, FluidsimSettings *domainSettings, Object *fsDomain,
        char *targetDir, char *targetFile)
{
	const char *suffixConfigTmp = FLUID_SUFFIX_CONFIG_TMP;

	/* prepare names... */
	const char *relbase = modifier_path_relbase(fsDomain);

	/* We do not accept empty paths, they can end in random places silently, see T51176. */
	if (domainSettings->surfdataPath[0] == '\0') {
		modifier_path_init(domainSettings->surfdataPath, sizeof(domainSettings->surfdataPath),
		                   OB_FLUIDSIM_SURF_DIR_DEFAULT);
		BKE_reportf(reports, RPT_WARNING, "Fluidsim: empty cache path, reset to default '%s'",
		            domainSettings->surfdataPath);
	}

	BLI_strncpy(targetDir, domainSettings->surfdataPath, FILE_MAXDIR);
	BLI_path_abs(targetDir, relbase);

	/* .tmp: don't overwrite/delete original file */
	BLI_join_dirfile(targetFile, FILE_MAX, targetDir, suffixConfigTmp);

	/* Ensure whole path exists and is wirtable. */
	const bool dir_exists = BLI_dir_create_recursive(targetDir);
	const bool is_writable = BLI_file_is_writable(targetFile);

	/* We change path to some presumably valid default value, but do not allow bake process to continue,
	 * this gives user chance to set manually another path. */
	if (!dir_exists || !is_writable) {
		modifier_path_init(domainSettings->surfdataPath, sizeof(domainSettings->surfdataPath),
		                   OB_FLUIDSIM_SURF_DIR_DEFAULT);

		if (!dir_exists) {
			BKE_reportf(reports, RPT_ERROR, "Fluidsim: could not create cache directory '%s', reset to default '%s'",
			            targetDir, domainSettings->surfdataPath);
		}
		else {
			BKE_reportf(reports, RPT_ERROR, "Fluidsim: cache directory '%s' is not writable, reset to default '%s'",
			            targetDir, domainSettings->surfdataPath);
		}

		BLI_strncpy(targetDir, domainSettings->surfdataPath, FILE_MAXDIR);
		BLI_path_abs(targetDir, relbase);

		/* .tmp: don't overwrite/delete original file */
		BLI_join_dirfile(targetFile, FILE_MAX, targetDir, suffixConfigTmp);

		/* Ensure whole path exists and is wirtable. */
		if (!BLI_dir_create_recursive(targetDir) || !BLI_file_is_writable(targetFile)) {
			BKE_reportf(reports, RPT_ERROR, "Fluidsim: could not use default cache directory '%s', "
			                                "please define a valid cache path manually", targetDir);
		}
		return false;
	}

	return true;
}

/* ******************************************************************************** */
/* ********************** write fluidsim config to file ************************* */
/* ******************************************************************************** */

typedef struct FluidBakeJob {
	/* from wmJob */
	void *owner;
	short *stop, *do_update;
	float *progress;
	int current_frame;
	elbeemSimulationSettings *settings;
} FluidBakeJob;

static void fluidbake_free(void *customdata)
{
	FluidBakeJob *fb= (FluidBakeJob *)customdata;
	MEM_freeN(fb);
}

/* called by fluidbake, only to check job 'stop' value */
static int fluidbake_breakjob(void *customdata)
{
	FluidBakeJob *fb= (FluidBakeJob *)customdata;

	if (fb->stop && *(fb->stop))
		return 1;
	
	/* this is not nice yet, need to make the jobs list template better 
	 * for identifying/acting upon various different jobs */
	/* but for now we'll reuse the render break... */
	return (G.is_break);
}

/* called by fluidbake, wmJob sends notifier */
static void fluidbake_updatejob(void *customdata, float progress)
{
	FluidBakeJob *fb= (FluidBakeJob *)customdata;
	
	*(fb->do_update) = true;
	*(fb->progress) = progress;
}

static void fluidbake_startjob(void *customdata, short *stop, short *do_update, float *progress)
{
	FluidBakeJob *fb= (FluidBakeJob *)customdata;
	
	fb->stop= stop;
	fb->do_update = do_update;
	fb->progress = progress;
	
	G.is_break = false;  /* XXX shared with render - replace with job 'stop' switch */
	
	elbeemSimulate();
	*do_update = true;
	*stop = 0;
}

static void fluidbake_endjob(void *customdata)
{
	FluidBakeJob *fb= (FluidBakeJob *)customdata;
	
	if (fb->settings) {
		MEM_freeN(fb->settings);
		fb->settings = NULL;
	}
}

static int runSimulationCallback(void *data, int status, int frame)
{
	FluidBakeJob *fb = (FluidBakeJob *)data;
	elbeemSimulationSettings *settings = fb->settings;
	
	if (status == FLUIDSIM_CBSTATUS_NEWFRAME) {
		fluidbake_updatejob(fb, frame / (float)settings->noOfFrames);
		//printf("elbeem blender cb s%d, f%d, domainid:%d noOfFrames: %d\n", status, frame, settings->domainId, settings->noOfFrames ); // DEBUG
	}
	
	if (fluidbake_breakjob(fb)) {
		return FLUIDSIM_CBRET_ABORT;
	}
	
	return FLUIDSIM_CBRET_CONTINUE;
}

static void fluidbake_free_data(FluidAnimChannels *channels, ListBase *fobjects, elbeemSimulationSettings *fsset, FluidBakeJob *fb)
{
	free_domain_channels(channels);
	MEM_freeN(channels);
	channels = NULL;

	free_all_fluidobject_channels(fobjects);
	BLI_freelistN(fobjects);
	MEM_freeN(fobjects);
	fobjects = NULL;
	
	if (fsset) {
		MEM_freeN(fsset);
		fsset = NULL;
	}
	
	if (fb) {
		MEM_freeN(fb);
		fb = NULL;
	}
}

/* copied from rna_fluidsim.c: fluidsim_find_lastframe() */
static void fluidsim_delete_until_lastframe(FluidsimSettings *fss, const char *relbase)
{
	char targetDir[FILE_MAX], targetFile[FILE_MAX];
	char targetDirVel[FILE_MAX], targetFileVel[FILE_MAX];
	char previewDir[FILE_MAX], previewFile[FILE_MAX];
	int curFrame = 1, exists = 0;

	BLI_join_dirfile(targetDir,    sizeof(targetDir),    fss->surfdataPath, OB_FLUIDSIM_SURF_FINAL_OBJ_FNAME);
	BLI_join_dirfile(targetDirVel, sizeof(targetDirVel), fss->surfdataPath, OB_FLUIDSIM_SURF_FINAL_VEL_FNAME);
	BLI_join_dirfile(previewDir,   sizeof(previewDir),   fss->surfdataPath, OB_FLUIDSIM_SURF_PREVIEW_OBJ_FNAME);

	BLI_path_abs(targetDir,    relbase);
	BLI_path_abs(targetDirVel, relbase);
	BLI_path_abs(previewDir,   relbase);

	do {
		BLI_strncpy(targetFile, targetDir, sizeof(targetFile));
		BLI_strncpy(targetFileVel, targetDirVel, sizeof(targetFileVel));
		BLI_strncpy(previewFile, previewDir, sizeof(previewFile));

		BLI_path_frame(targetFile, curFrame, 0);
		BLI_path_frame(targetFileVel, curFrame, 0);
		BLI_path_frame(previewFile, curFrame, 0);

		curFrame++;

		if ((exists = BLI_exists(targetFile))) {
			BLI_delete(targetFile, false, false);
			BLI_delete(targetFileVel, false, false);
			BLI_delete(previewFile, false, false);
		}
	} while (exists);

	return;
}

static int fluidsimBake(bContext *C, ReportList *reports, Object *fsDomain, short do_job)
{
	Scene *scene = CTX_data_scene(C);
	ViewLayer *view_layer = CTX_data_view_layer(C);
	Depsgraph *depsgraph = CTX_data_depsgraph(C);
	int i;
	FluidsimSettings *domainSettings;

	char debugStrBuffer[256];
	
	int gridlevels = 0;
	const char *relbase= modifier_path_relbase(fsDomain);
	const char *strEnvName = "BLENDER_ELBEEMDEBUG"; // from blendercall.cpp
	const char *suffixConfigTmp = FLUID_SUFFIX_CONFIG_TMP;
	const char *suffixSurface = FLUID_SUFFIX_SURFACE;

	char targetDir[FILE_MAX];  // store & modify output settings
	char targetFile[FILE_MAX]; // temp. store filename from targetDir for access

	float domainMat[4][4];
	float invDomMat[4][4];

	int noFrames;
	int origFrame = scene->r.cfra;
	
	FluidAnimChannels *channels = MEM_callocN(sizeof(FluidAnimChannels), "fluid domain animation channels");
	ListBase *fobjects = MEM_callocN(sizeof(ListBase), "fluid objects");
	FluidsimModifierData *fluidmd = NULL;
	Mesh *mesh = NULL;

	FluidBakeJob *fb;
	elbeemSimulationSettings *fsset= MEM_callocN(sizeof(elbeemSimulationSettings), "Fluid sim settings");

	fb= MEM_callocN(sizeof(FluidBakeJob), "fluid bake job");
	
	if (getenv(strEnvName)) {
		int dlevel = atoi(getenv(strEnvName));
		elbeemSetDebugLevel(dlevel);
		BLI_snprintf(debugStrBuffer, sizeof(debugStrBuffer), "fluidsimBake::msg: Debug messages activated due to envvar '%s'\n", strEnvName);
		elbeemDebugOut(debugStrBuffer);
	}
	
	/* make sure it corresponds to startFrame setting (old: noFrames = scene->r.efra - scene->r.sfra +1) */;
	noFrames = scene->r.efra - 0;
	if (noFrames<=0) {
		BKE_report(reports, RPT_ERROR, "No frames to export (check your animation range settings)");
		fluidbake_free_data(channels, fobjects, fsset, fb);
		return 0;
	}
	
	/* check scene for sane object/modifier settings */
	if (!fluid_validate_scene(reports, view_layer, fsDomain)) {
		fluidbake_free_data(channels, fobjects, fsset, fb);
		return 0;
	}
	
	/* these both have to be valid, otherwise we wouldn't be here */
	fluidmd = (FluidsimModifierData *)modifiers_findByType(fsDomain, eModifierType_Fluidsim);
	domainSettings = fluidmd->fss;
	mesh = fsDomain->data;
	
	domainSettings->bakeStart = 1;
	domainSettings->bakeEnd = scene->r.efra;
	
	// calculate bounding box
	fluid_get_bb(mesh->mvert, mesh->totvert, fsDomain->obmat, domainSettings->bbStart, domainSettings->bbSize);
	
	// reset last valid frame
	domainSettings->lastgoodframe = -1;

	/* delete old baked files */
	fluidsim_delete_until_lastframe(domainSettings, relbase);
	
	/* rough check of settings... */
	if (domainSettings->previewresxyz > domainSettings->resolutionxyz) {
		BLI_snprintf(debugStrBuffer, sizeof(debugStrBuffer), "fluidsimBake::warning - Preview (%d) >= Resolution (%d)... setting equal.\n", domainSettings->previewresxyz,  domainSettings->resolutionxyz);
		elbeemDebugOut(debugStrBuffer);
		domainSettings->previewresxyz = domainSettings->resolutionxyz;
	}
	// set adaptive coarsening according to resolutionxyz
	// this should do as an approximation, with in/outflow
	// doing this more accurate would be overkill
	// perhaps add manual setting?
	if (domainSettings->maxRefine <0) {
		if (domainSettings->resolutionxyz>128) {
			gridlevels = 2;
		}
		else if (domainSettings->resolutionxyz > 64) {
			gridlevels = 1;
		}
		else {
			gridlevels = 0;
		}
	}
	else {
		gridlevels = domainSettings->maxRefine;
	}
	BLI_snprintf(debugStrBuffer, sizeof(debugStrBuffer), "fluidsimBake::msg: Baking %s, refine: %d\n", fsDomain->id.name, gridlevels);
	elbeemDebugOut(debugStrBuffer);
	
	
	
	/* ******** prepare output file paths ******** */
	if (!fluid_init_filepaths(reports, domainSettings, fsDomain, targetDir, targetFile)) {
		fluidbake_free_data(channels, fobjects, fsset, fb);
		return false;
	}

	channels->length = scene->r.efra; // DG TODO: why using endframe and not "noFrames" here? .. because "noFrames" is buggy too? (not using sfra)
	channels->aniFrameTime = (double)((double)domainSettings->animEnd - (double)domainSettings->animStart) / (double)noFrames;
	
	/* ******** initialize and allocate animation channels ******** */
	fluid_init_all_channels(C, fsDomain, domainSettings, channels, fobjects);

	/* reset to original current frame */
	scene->r.cfra = origFrame;
	ED_update_for_newframe(CTX_data_main(C), depsgraph);
		
	/* ******** init domain object's matrix ******** */
	copy_m4_m4(domainMat, fsDomain->obmat);
	if (!invert_m4_m4(invDomMat, domainMat)) {
		BLI_snprintf(debugStrBuffer, sizeof(debugStrBuffer), "fluidsimBake::error - Invalid obj matrix?\n");
		elbeemDebugOut(debugStrBuffer);
		BKE_report(reports, RPT_ERROR, "Invalid object matrix"); 

		fluidbake_free_data(channels, fobjects, fsset, fb);
		return 0;
	}

	/* ********  start writing / exporting ******** */
	// use .tmp, don't overwrite/delete original file
	BLI_join_dirfile(targetFile, sizeof(targetFile), targetDir, suffixConfigTmp);

	/* ******** export domain to elbeem ******** */
	elbeemResetSettings(fsset);
	fsset->version = 1;
	fsset->threads = (domainSettings->threads == 0) ? BKE_scene_num_threads(scene) : domainSettings->threads;
	// setup global settings
	copy_v3_v3(fsset->geoStart, domainSettings->bbStart);
	copy_v3_v3(fsset->geoSize, domainSettings->bbSize);
	
	// simulate with 50^3
	fsset->resolutionxyz = (int)domainSettings->resolutionxyz;
	fsset->previewresxyz = (int)domainSettings->previewresxyz;

	fsset->realsize = get_fluid_size_m(scene, fsDomain, domainSettings);
	fsset->viscosity = get_fluid_viscosity(domainSettings);
	get_fluid_gravity(fsset->gravity, scene, domainSettings);

	// simulate 5 frames, each 0.03 seconds, output to ./apitest_XXX.bobj.gz
	fsset->animStart = domainSettings->animStart;
	fsset->aniFrameTime = channels->aniFrameTime;
	fsset->noOfFrames = noFrames; // is otherwise subtracted in parser

	BLI_join_dirfile(targetFile, sizeof(targetFile), targetDir, suffixSurface);

	// defaults for compressibility and adaptive grids
	fsset->gstar = domainSettings->gstar;
	fsset->maxRefine = domainSettings->maxRefine; // check <-> gridlevels
	fsset->generateParticles = domainSettings->generateParticles; 
	fsset->numTracerParticles = domainSettings->generateTracers; 
	fsset->surfaceSmoothing = domainSettings->surfaceSmoothing; 
	fsset->surfaceSubdivs = domainSettings->surfaceSubdivs; 
	fsset->farFieldSize = domainSettings->farFieldSize; 
	BLI_strncpy(fsset->outputPath, targetFile, sizeof(fsset->outputPath));

	// domain channels
	fsset->channelSizeFrameTime = 
	fsset->channelSizeViscosity = 
	fsset->channelSizeGravity = channels->length;
	fsset->channelFrameTime = channels->DomainTime;
	fsset->channelViscosity = channels->DomainViscosity;
	fsset->channelGravity = channels->DomainGravity;
	
	fsset->runsimCallback = &runSimulationCallback;
	fsset->runsimUserData = fb;

	if (domainSettings->typeFlags & OB_FSBND_NOSLIP)		fsset->domainobsType = FLUIDSIM_OBSTACLE_NOSLIP;
	else if (domainSettings->typeFlags&OB_FSBND_PARTSLIP)	fsset->domainobsType = FLUIDSIM_OBSTACLE_PARTSLIP;
	else if (domainSettings->typeFlags&OB_FSBND_FREESLIP)	fsset->domainobsType = FLUIDSIM_OBSTACLE_FREESLIP;
	fsset->domainobsPartslip = domainSettings->partSlipValue;

	/* use domainobsType also for surface generation flag (bit: >=64) */
	if (domainSettings->typeFlags & OB_FSSG_NOOBS)
		fsset->mFsSurfGenSetting = FLUIDSIM_FSSG_NOOBS;
	else
		fsset->mFsSurfGenSetting = 0; // "normal" mode

	fsset->generateVertexVectors = (domainSettings->domainNovecgen==0);

	// init blender domain transform matrix
	{ int j;
	for (i=0; i<4; i++) {
		for (j=0; j<4; j++) {
			fsset->surfaceTrafo[i*4+j] = invDomMat[j][i];
		}
	} }

	/* ******** init solver with settings ******** */
	elbeemInit();
	elbeemAddDomain(fsset);
	
	/* ******** export all fluid objects to elbeem ******** */
	export_fluid_objects(C, fobjects, scene, channels->length);
	
	/* custom data for fluid bake job */
	fb->settings = fsset;
	
	if (do_job) {
		wmJob *wm_job = WM_jobs_get(CTX_wm_manager(C), CTX_wm_window(C), scene, "Fluid Simulation",
		                            WM_JOB_PROGRESS, WM_JOB_TYPE_OBJECT_SIM_FLUID);

		/* setup job */
		WM_jobs_customdata_set(wm_job, fb, fluidbake_free);
		WM_jobs_timer(wm_job, 0.1, NC_SCENE|ND_FRAME, NC_SCENE|ND_FRAME);
		WM_jobs_callbacks(wm_job, fluidbake_startjob, NULL, NULL, fluidbake_endjob);

		WM_jobs_start(CTX_wm_manager(C), wm_job);
	}
	else {
		short dummy_stop = 0, dummy_do_update = 0;
		float dummy_progress = 0.0f;

		/* blocking, use with exec() */
		fluidbake_startjob((void *)fb, &dummy_stop, &dummy_do_update, &dummy_progress);
		fluidbake_endjob((void *)fb);
		fluidbake_free((void *)fb);
	}

	/* ******** free stored animation data ******** */
	fluidbake_free_data(channels, fobjects, NULL, NULL);

	// elbeemFree();
	return 1;
}

static void UNUSED_FUNCTION(fluidsimFreeBake)(Object *UNUSED(ob))
{
	/* not implemented yet */
}

#else /* WITH_MOD_FLUID */

/* only compile dummy functions */
static int fluidsimBake(bContext *UNUSED(C), ReportList *UNUSED(reports), Object *UNUSED(ob), short UNUSED(do_job))
{
	return 0;
}

#endif /* WITH_MOD_FLUID */

/***************************** Operators ******************************/

static int fluid_bake_invoke(bContext *C, wmOperator *op, const wmEvent *UNUSED(event))
{
	/* only one bake job at a time */
	if (WM_jobs_test(CTX_wm_manager(C), CTX_data_scene(C), WM_JOB_TYPE_OBJECT_SIM_FLUID))
		return OPERATOR_CANCELLED;

	if (!fluidsimBake(C, op->reports, CTX_data_active_object(C), true))
		return OPERATOR_CANCELLED;

	return OPERATOR_FINISHED;
}

static int fluid_bake_exec(bContext *C, wmOperator *op)
{
	if (!fluidsimBake(C, op->reports, CTX_data_active_object(C), false))
		return OPERATOR_CANCELLED;

	return OPERATOR_FINISHED;
}

void FLUID_OT_bake(wmOperatorType *ot)
{
	/* identifiers */
	ot->name = "Fluid Simulation Bake";
	ot->description = "Bake fluid simulation";
	ot->idname = "FLUID_OT_bake";
	
	/* api callbacks */
	ot->invoke = fluid_bake_invoke;
	ot->exec = fluid_bake_exec;
	ot->poll = ED_operator_object_active_editable;
}