Updated Melting

After spring break and diving back into this project, I've been working on stabilizing and trying to reduce numerical issues/errors with the melting simulation as it currently exists.

(The paper I will refer to below is here)

http://www.matmidia.mat.puc-rio.br/tomlew/pdfs/plastic_cad.pdf

It seemed like my volume loss during melting was due to a couple of errors in the way I implemented the paper, as well as a poor selection of constants. The problem with melting is that we add several new constants to control viscosity of the fluid, which need some thorough tweaking to play nicely with the rest of the simulation. 

As far as other errors, I decided to delve back into my computational linear textbook to try to figure out if I was computing some of the tensor calculations incorrectly. All I will say here is that a tensor dot product with a vector (T . m) is not the same as a matrix M x v (see http://www.polymerprocessing.com/notes/root92a.pdf).  After this I was able to manually set the viscosity and have a lot more control of the melting situation control. The paper describes a "jump" number which also affects the viscosity according to the local shear of particles, but I have not been able to extract good viscosity values from the function they use and so I am unable to simulate solids using that jump number. 

While I have tied my thermal diffusion of temperature function successfully to this jump number which is producing reasonable values, the jump number to viscosity function which is also dependent on deformation tensor magnitude (essentially "shear") is producing viscosity values which are way too low. Because the function is non-linear I do not think I can just scale it up so I will need to either figure out if there is an error or find a new method of viscosity adjustment. Another option would be to ditch the "shear" value altogether and have visc. vary linearly with temperature and the "jump" value. I'm not sure what results this will produce but it is worth a try. Any thoughts on what I should do in this case?

If the simulation is too viscous, the simulation is unstable at the surface boundary so we need to decrease the timestep (the paper claims it solves this by approaching the pressure solve by taking the pressure derivative then updating it using Leap Frog integration, but after testing this method out I actually found it de-stablized my simulation even more and led to more unwanted compression during certain simulations). However I did find a reasonably small timestep which doesn't explode and keeps rigid bodies stable. 

Here is a recap for this week:

(a) Stablized updating stress/deformation tensors through linear algebra magic.

(b) Tested out different methods of updating pressure forces.

• Decided on first-order density update. SPH pressure force approximation using density.

(c) Implemented XSPH particle velocity correction

• Adjusts particle velocities with an average velocity from velocities of neighboring particles to avoid  stable clumps of particles.

(d)  New and Improved Localized Artificial Viscosity!

• To avoid numerical instabilities in the velocity vector field, we adjust the acceleration using a different method which takes into effect local velocities and position differences when adjusting acceleration for art. visc. 

Rigid Body cube being pretty stable after about 1500 iterations (its a guess but it was running for like 5 minutes at around 4 fps so I'm going to guess it was that much...)

Same simulation run at alower viscosity.

Same simulation run at much lower viscosity. Notice the green which shows that area has a higher velocity in contrast with the other simulations.