← Back to team overview

dolfin team mailing list archive

Re: vector elements/single dof imposed

 

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
((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



--
Vuk Milisic
Chercheur CNRS
UMR 5224
Laboratoire Jean Kuntzman
Unversité Joseph Fourier
Grenoble
http://ljk.imag.fr/membres/Vuk.Milisic/


Follow ups

References