dolfin team mailing list archive

[Anders Logg <logg@xxxxxxxxx>]

On Thu, Jun 26, 2008 at 11:44:58PM +0200, Vuk.Milisic@xxxxxxx wrote:
> ok I followed your advice and set in the ffc part:
> 
> P2 = FiniteElement("Lagrange", "triangle", 2)
> P1 = FiniteElement("Lagrange", "triangle", 1)
> Vh = P2 + P2 +P2
> Th = (P2 + P2 +P2) + P1

You should not need to define the vector-valued space as a mixed
space, just do

  P2 = VectorElement("Lagrange", "triangle", 2)
       ~~~~~~

Writing

  Vh = P2 + P2 + P2

does not do what you might expect. What happens is that first
P2 + P2 yields a mixed element (with two components), then
(P2 + P2) + P2 yields a nested mixed element with two subsystems,
the first being P2 + P2 and the second being P2.

The correct way to do this (if you don't want to use VectorElement
which is the simplex way) is

  Vh = MixedElement([P2, P2, P2])

-- 
Anders


> ((v0, v1, v2), q) = TestFunctions(Th)
> ((u0, u1, u2), p) = TrialFunctions(Th)
> ((phi0, phi1, phi2),pkmu) = Functions(Th)
> 
> 
> then I wrote bc on single components
> and set
> 
> 
>    SubSystem velocity0(0,0);
>    SubSystem velocity1(0,1);
>    SubSystem velocity2(0,2);
>    SubSystem pressure (1);
>    // Set up boundary condition at left end
>    zero c(mesh);
>    Left left;
>    DirichletBC bcl(c, mesh, left,velocity0);
> 
>    // Set up boundary condition at right end
>    Pull0 r(mesh);
>    Right right;
>    DirichletBC bcr0(r, mesh, right,velocity0);
>    DirichletBC bcr1(c, mesh, right,velocity1);
> 
>    Bottom bottom;
>    DirichletBC bcb1(c, mesh,bottom,velocity1 );
>    DirichletBC bcb2(c, mesh,bottom,velocity2 );
> 
>    // Set up boundary conditions
>    Array<BoundaryCondition*> bcs;
>    bcs.push_back(&bcl);
>    bcs.push_back(&bcr0);
>    bcs.push_back(&bcr1);
>    bcs.push_back(&bcb1);
>    bcs.push_back(&bcb2);
> 
> 
> 
> The compile is ok, the runtime gives:
> 
> Creating sub domain markers for boundary condition.
> Creating sub domain markers for boundary condition.
> Creating sub domain markers for boundary condition.
> Creating sub domain markers for boundary condition.
> Creating sub domain markers for boundary condition.
> Creating linear PDE with 5 boundary condition(s).
> Solving linear PDE.
>    | Assembling matrix over cells                                    |
>    |===============================|---------------------------------| 47.9%
>    | Assembling matrix over cells                                    |
>    |==============================================================|--| 96.1%
>    | Assembling matrix over cells                                    |
>    |=================================================================| 100.0%
>    Assembling vector over cells (finished).
>    Computing Dirichlet boundary values (finished).
>    Applying boundary conditions to linear system.
>    Computing Dirichlet boundary values (finished).
>    Applying boundary conditions to linear system.
>    Computing Dirichlet boundary values (finished).
>    Applying boundary conditions to linear system.
>    Computing Dirichlet boundary values (finished).
>    Applying boundary conditions to linear system.
> terminate called after throwing an instance of 'std::runtime_error'
>    what():  *** Error: Unable to extract sub system 2 (only 2 sub  
> systems defined).
> Abandon
> 
> 
> It seems that he does not reach the dofs associated to u2 v2 phi2
> when it should, is it a bug or I am missing something about the
> boundary condition ?
> 
> 
> 
> Cheers,
> 
> 
> VM
> Quoting "Garth N. Wells" <gnw20@xxxxxxxxx>:
> 
> > Vuk Milisic wrote:
> >> If I want to fix a boundary condition on a single component of a   
> >> vector element
> >> u[0] and not touch the other components how should I define the heritated
> >> class of Function ?
> >>
> >
> > Look at the Stokes demo,
> >
> > demo/pde/stokes/taylor-hood/cpp
> >
> > Garth
> >
> >
> >>
> >> if I do
> >> /class Slide : public Function
> >>  {
> >>  public:
> >>
> >>    Slide(Mesh& mesh) : Function(mesh) {}
> >>
> >>    void eval(real* values, const real* x) const
> >>    {
> >>      //      values[0] = 0.0;
> >>      values[1] = 0.0;
> >>      values[2] = 0.0;
> >>    }
> >>
> >>  };/
> >> It fixes also the first component. How should I proceed ?
> >>
> >> Thanks,
> >>
> >>
> >> VM
> >>
> >> -- 
> >>
> >>
> >>
> >> -----------------------------------
> >> Vuk Milisic
> >> Chargé de Recherche
> >> LJK-IMAG UMR 5523
> >> 51, rue des Mathematiques -  B. P. 53 38041 Grenoble Cedex 9    
> >> France Office   n° 66
> >> Tel:  +33 4 76 63 57 38
> >> Fax:  +33 4 76 63 12 63
> >> http://ljk.imag.fr/membres/Vuk.Milisic/
> >> -----------------------------------
> >>
> >>
> >> ------------------------------------------------------------------------
> >>
> >> _______________________________________________
> >> DOLFIN-dev mailing list
> >> DOLFIN-dev@xxxxxxxxxx
> >> http://www.fenics.org/mailman/listinfo/dolfin-dev
> 
> 
>

Attachment: signature.asc
Description: Digital signature