dolfin team mailing list archive

["Johan Hoffman" <jhoffman@xxxxxxxxxx>]

> On Tue, May 30, 2006 at 10:07:47AM +0000, Alexander Jarosch wrote:
>> Hi everybody.
>>
>> There was a short discussion about Stokes flow on complex geometries. I
>> would like to report that a ffc form like the one below works well on
>> complex geometries, I still only have a factor of 1.2 offset in
>> magnitude if I compare it with another fem model.  I would like to get
>> some feedback about this form, whether it is okay or not from the
>> developers point of view.
>>
>> much appreciated,
>>
>> Alex
>>
>> ffc form -------------------
>>
>> scalar = FiniteElement("Lagrange", "triangle", 1)
>> vector = FiniteElement("Vector Lagrange", "triangle", 1)
>> system = vector + scalar
>>
>> (v, q) = TestFunctions(system)
>> (u, p) = TrialFunctions(system)
>>
>> f = Function(vector)
>> h = Function(scalar)
>> nu = Function(scalar)
>>
>> beta  = 0.2
>> delta = beta*h*h
>>
>> a = (nu*dot(grad(v), grad(u)) - div(v)*p + nu*q*div(u) +
>> delta*dot(grad(q), grad(p)))*dx
>> L = dot(v + mult(delta, grad(q)), f)*dx
>
> I think it looks a little strange. The difference from what you had
> before seems to be that you have included the viscosity in the
> incompressibility term:
>
>     nu*q*div(u)
>
> Is this what makes the difference? If so, I don't know why this
> helps. It's certainly ok to put it there since you're solving for
> div(u) = 0 anyway, so nu*div(u) = 0 is also correct.
>
> In some way, it corresponds to a smaller penalty for compressibility
> when the viscosity is small which might be reasonable.
>
> /Anders

I think this term is somewhat strange, and I haven't seen it before. I do
not see why the viscosity should multiply the incompressibility.

If the problem is connected to the stability of the incompressibility
constraint, I propose to instead add to the left hand side of the form a
stabilization term for the incompressibility, of the form:

delta2*div(v)*div(u)

with delta2 = C2*h*h, and C2 of unit order.

Also, what is the value of your parameter beta in the definition of your
delta? It might be worthhile to try to experiment with higher values of
beta.

/Johan



>
> _______________________________________________
> DOLFIN-dev mailing list
> DOLFIN-dev@xxxxxxxxxx
> http://www.fenics.org/cgi-bin/mailman/listinfo/dolfin-dev
>