Method init

SPlisHSPlasH/PCISPH/SimulationDataPCISPH.cpp:19–108

}

void SimulationDataPCISPH::init()
{
	Simulation *sim = Simulation::getCurrent();
	const unsigned int nModels = sim->numberOfFluidModels();

	m_predX.resize(nModels);
	m_predV.resize(nModels);
	m_densityAdv.resize(nModels);
	m_pressure.resize(nModels);
	m_pressureAccel.resize(nModels);
	m_pcisph_factor.resize(nModels);
	for (unsigned int i = 0; i < nModels; i++)
	{
		FluidModel *fm = sim->getFluidModel(i);
		m_predX[i].resize(fm->numParticles(), Vector3r::Zero());
		m_predV[i].resize(fm->numParticles(), Vector3r::Zero());
		m_densityAdv[i].resize(fm->numParticles(), 0.0);
		m_pressure[i].resize(fm->numParticles(), 0.0);
		m_pressureAccel[i].resize(fm->numParticles(), Vector3r::Zero());
	}

	LOG_INFO << "Initialize PCISPH scaling factor";
	for (unsigned int fluidModelIndex = 0; fluidModelIndex < nModels; fluidModelIndex++)
	{
		FluidModel *model = sim->getFluidModel(fluidModelIndex);
		m_pcisph_factor[fluidModelIndex] = 0.0;

		// Find prototype particle
		// => particle with max. fluid neighbors
		const Real density0 = model->getDensity0();
 
 		Vector3r sumGradW = Vector3r::Zero();
 		Real sumGradW2 = 0.0;
 		const Real supportRadius = sim->getSupportRadius();
 		const Real particleRadius = sim->getParticleRadius();
 		const Real diam = static_cast<Real>(2.0) * particleRadius;
 		const Vector3r xi(0,0,0);
 
 		// use a regular sampling around (0,0,0)
		if (sim->is2DSimulation())
		{
			Vector3r xj = { -supportRadius, -supportRadius, 0.0 };
			while (xj[0] <= supportRadius)
			{
				while (xj[1] <= supportRadius)
				{
					// check if xj is in the support of xi
					if ((xi - xj).squaredNorm() < supportRadius*supportRadius)
					{
						const Vector3r gradW = sim->gradW(xi - xj);
						sumGradW += gradW;
						sumGradW2 += gradW.squaredNorm();
					}
					xj[1] += diam;
				}
				xj[0] += diam;
				xj[1] = -supportRadius;
			}