dolfin team mailing list archive

["Garth N. Wells" <g.n.wells@xxxxxxxxxx>]

On Wed, 2006-02-22 at 10:45 +0000, Alexander Jarosch wrote:
> Hi,
> 
> I have a rather simple question I guess. I tried to figure out how to 
> give the stokes solver a viscosity for each vertex and started of with 
> the stokes solver with a constant viscosity nu:
> 
> modified from Anders file:
> ---------------
> P1 = FiniteElement("Lagrange", "triangle", 1)
> P2 = FiniteElement("Vector Lagrange", "triangle", 2)
> TH = P2 + P1
> 
> (v, q) = BasisFunctions(TH)
> (u, p) = BasisFunctions(TH)
> nu     = Constant()    #viscosity
> 
> 
> f = Function(P2)
> 
> a = (dot(nu*grad(u), grad(v)) - p*div(v) + div(u)*q)*dx
> L = dot(f, v)*dx
> ---------------
> 
> now that works fine and now i have the viscosity values for each vertex 
> in a vector as values and in a function like:
> 
> Function viscfunc(viscosity, mesh, element);
> 
> Now i tried to get a form file which takes these values like:
> ---------------
> P1 = FiniteElement("Lagrange", "triangle", 1)
> P2 = FiniteElement("Vector Lagrange", "triangle", 2)
> element1 = FiniteElement("Discontinuous vector Lagrange", "triangle", 0, 1)
> 

element1 is a constant, so you won't be able to set the viscosity at
vertexes, just on elements. Also, you can simplify it by writing
 
element1 = FiniteElement("Discontinuous Lagrange", "triangle", 0)

> 
> 
> 
> TH = P2 + P1
> 
> (v, q) = BasisFunctions(TH)
> (u, p) = BasisFunctions(TH)
> nu     = Function(element1)    #viscosity
> 
> f = Function(P2)
> 
> a = (dot(dot(nu, grad(u)), grad(v)) - p*div(v) + div(u)*q)*dx
> L = dot(f, v)*dx
> ---------------
> 
> well and it does not work, as well as the approach like:
> 
> ---------------
> P1 = FiniteElement("Lagrange", "triangle", 1)
> P2 = FiniteElement("Vector Lagrange", "triangle", 2)
> constant_scalar = FiniteElement("Discontinuous Lagrange", "triangle", 1)
> 
> 

"constant_scalar" is not constant. It's linear on elements, but because
it's discontinuous it's value is not unique at each vertex in the mesh.
Try 

constant_scalar = FiniteElement("Discontinuous Lagrange", "triangle", 0)

or

linear_scalar = FiniteElement("Lagrange", "triangle", 1)


> 
> TH = P2 + P1
> 
> (v, q) = BasisFunctions(TH)
> (u, p) = BasisFunctions(TH)
> nu     = Function(constant_scalar)    #viscosity
> #nu     = Constant()
> 
> f = Function(P2)
> 
> a = (dot(nu*grad(u), grad(v)) - p*div(v) + div(u)*q)*dx
> L = dot(f, v)*dx
> ---------------
> 
> which actually takes my viscosity vector and starts to solve but does 
> not converge.

Make sure that you initialise Function correctly in DOLFIN. Try with a
constant function first,

Function nu = 1.0;

Garth

> 
> 
> 
> 
> 
> might be a simple mistake in there, I appreciate any help or hints.
> 
> Alex
> 
> _______________________________________________
> DOLFIN-dev mailing list
> DOLFIN-dev@xxxxxxxxxx
> http://www.fenics.org/cgi-bin/mailman/listinfo/dolfin-dev
-- 
Dr. Garth N. Wells
Faculty of Civil Engineering and Geosciences
Delft University of Technology
Stevinweg 1
2628 CN Delft
The Netherlands

tel.     +31 15 278 7922
fax.     +31 15 278 6383
e-mail   g.n.wells@xxxxxxxxxx
url      http://www.mechanics.citg.tudelft.nl/~garth