Discussion Forum

I am trying to input a viscous model for a squeeze flow simulation where it depends not only on the second invariant but also in the third invariant

I noted that comsol inputs the viscosity as a function of the shear rate which is a function of the gradients(ux,uy,vx,uz,wx,vy,vz,wz,wy) in fact the shear rate is completely a scalar value of the second invariant as normally is simulated in fluid dynamics


in my case my model not only depends on the shear rate depends on the second and the third invariant is a viscosity function like this
η(II,III)=η_0 (III/(6(II/6)^(3/2) ))+m(√((((II)-6(III/6)^(2/3) ))/2))^(n-1)


so when I input in the user defined viscosity a function of the third and second invariant as a functions of the gradients(ux,uy,vx,uz,wx,vy,vz,wz,wy), and I try to solve the model I encounter with several erros that I dont know if are part of my function input, part of the solver or part the model

My first try was to input the model in 2D to make it simplified and with a model just as a function of the firts term without taking in account the exponential term there i found the error: "Failed to find a solution.Maximum number of Newton iterations reached.Returned solution is not converged."

the viscosity function was this:
100000[Pa*s]*(((uy+vx+2*vx)*((uy+vx)^2+(2*vx)^2 )+(uy+vx)^2*(2*ux+2*vx)+2*ux*((2*ux)^2+(uy+vx)^2 )+ eps)/(6*(((2*ux)^2+2*(uy+vx)^2+(2*vx)^2+eps)/6)^(3/2) ))+eps

which is η(II,III)=η_0 (III/(6(II/6)^(3/2) ))

then when i tried to input it in a 3D model with the same viscosity function to be like:
10000[Pa*s]*(((2*uy + 2*vx)*(uz + wx)*(vz + wy)- 2*ux*(vz + wy)^2- 2*wz*(uy + vx)^2- 2*vx*(uz + wx)^2+ 8*ux*vx*wz+eps)/(6*(((2*ux)^2+2*(vx+uy)^2+2*(wy+vz)^2+2*(uz+wx)^2+(2*vy)^2+(2*wz)^2+eps)/6)^(3/2) ))+eps

and I found the next kind of error:
Detail: NaN or Inf found when solving linear system using SOR.

evenmore when I tried to substituing the complete viscous model like this:

η(II,III)=η_0 ((5((vx+uy)*(wy+vz)*(uz+wx))+2(wy+vz)^2*(2wz+2vy)+(vx+uy+2vy)*((vx+uy)^2+(2vy)^2+(wy+vz)^2 )+2(uz+wx)^2*(2ux+2wz)+(uy+vx)^2*(2ux+2vy)+(2ux)*((2ux)^2+(uy+vx)^2+(uz+wx)^2 )+(2wz)*((2wz)^2+(uz+wx)^2+(vz+wy)^2 ))/(6(((2ux)^2+2(vx+uy)^2+2(wy+vz)^2+2(uz+wx)^2+(2vy)^2+(2wz)^2)/6)^(3/2) ))+m(√(((((2ux)^2+2(vx+uy)^2+2(wy+vz)^2+2(uz+wx)^2+(2vy)^2+(2wz)^2 )-6((5((vx+uy)*(wy+vz)*(uz+wx))+2(wy+vz)^2*(2wz+2vy)+(vx+uy+2vy)*((vx+uy)^2+(2vy)^2+(wy+vz)^2 )+2(uz+wx)^2*(2ux+2wz)+(uy+vx)^2*(2ux+2vy)+(2ux)*((2ux)^2+(uy+vx)^2+(uz+wx)^2 )+(2wz)*((2wz)^2+(uz+wx)^2+(vz+wy)^2 ))/6)^(2/3) ))/2))^(n-1)

I got an error like:"Attempt to evaluate real square root of negative number. Function: sqrt"

an interesting thing is that for the same model, with the same boundary conditions I just found a solution when the model only depends on the 2 invariant with a viscosity input like that:

η(II,III)=m*(sqrt(((((2*ux)^2+2*(vx+uy)^2+2*(wy+vz)^2+2*(uz+wx)^2+(2*vy)^2+(2*wz)^2+eps) ))/2))^(n-1)

So I am wondering, where my errors come from and how I can get a result using the third invariant or what kind of recommendations I can receive to input this viscosity model

I attached the files in case they are needed

Thanks Roberto