Re: [Getfem-users] vectorial problem with interpolated fem

[Andriy Andreykiv - 3ME]

Dear Yves, 

Some time ago you asked me (see below)  to check if the 
getfem::interpolated_fem works, namely the basis functions and the gradients. 
Back then I said the basis functions were working fine and now I can confirm 
that the gradients seem to work correctly as well.

A question on another mesh_fem object, if I may. I'm developping a contact 
algorithm that uses Lagrange multipliers to enforce contact constraints.
First I tried to build a mesh on the boundary of the master body for the 
Lagrange multiplier space and remove the elements (lines in 2D or faces in 3D) 
located ouside contact area.
The idea is OK, but when I remove elements of this LM mesh, I spoil the 
numbering of the DOF's and the element numbers as well. So, for instance, If I 
want to remove elements 1 and 2 from 
the mesh, I would remove element 1, but the mesh would renumber and the element 
2 is already not the one that I want to remove.

The second idea was to use this new class that you developed probably 
specifically for such cases: getfem::partial_mesh_fem.
But it still has this limitation that it doesn't allow vectorial FEMs (I got an 
error about it when I ran my program). 
And I do need vectorial FEM's to calculate, say, contact forces in the nodes 
and the corresponding stiffness matrix.

Does Getfem allow some easy solution to this problem, or I'll have to reprogram 
my code so it doesn't use vectorial fems?

Thank you very much in advance,
                                    Andriy



 

-----Original Message-----
From: Yves Renard [mailto:address@hidden
Sent: donderdag 29 januari 2009 10:04
To: Andriy Andreykiv - 3ME
Subject: Re: vectorial problem with interpolated fem


Dear Andriy,

I checked my changes and I found an error. The new file is attached. I think it 
is now ok at least for the base functions. If you could carrefully test for the 
gradients it would be nice. If it is ok, I commit this change for the future 
versions of Getfem.

Yves.




On Tuesday 27 January 2009 21:21, you wrote:
> Dear Yves,
>
> Thank you. Sorry for the second e-mail yesterday, I had problems with 
> e-mail server and couldn't see you that wrote me the first time.
>
> Your suggestion worked and I passed the point of interpolate fem 
> creation, however, now it gives another error. Right after the last 
> line in the code below,
>
> //mf_interpolate - mesh_fem from interpolation of 2D on a 1D mesh, qdim=1
> //mf_U_line      - mesh_fem for that 1D line with qdim=2
> //mf_U_interpolate - the same as mf_interpolate but with qdim=2 (our vector
> interpolated fem) //mf_lambda        - mesh_fem for that line with qdim=1
>
>     LAMBDA.resize(mf_lambda.nb_dof());
>    
> 
>gmm::copy(gmm::sub_vector(LAMBDA_FULL,gmm::unsorted_sub_index(index_imp
>l_la mbda)),LAMBDA); plain_vector LS_small(mf_interpolate.nb_dof());
>     plain_vector RHS_L_ui(mf_U_line.nb_dof());
>     plain_vector RHS_L_ub(mf_U_interpolate.nb_dof());
>     plain_vector RHS_L_L(mf_lambda.nb_dof());
>     gmm::copy(gmm::sub_vector(ls.values(),
>                 gmm::unsorted_sub_index(index_intpl_ls)),LS_small);
>     getfem::generic_assembly assem_boundary;
>     assem_boundary.set(    "F=data$1(#3);"
>                          "L=data$2(#1);"
>                          
> "V$1(#2)+=comp(Base(#1).Grad(#3).vBase(#2))(i,j,k,:,k).L(i).F(j);"
>                         
> "V$2(#4)+=(-1.0)*comp(Base(#1).Grad(#3).vBase(#4))(i,j,k,:,k).L(i).F(j);"
> "V$3(#1)+=comp(Base(#1).Base(#3)) (:,i).F(i);"
>                            );
>     assem_boundary.push_mi(mim_line);        //mim on the line
>     assem_boundary.push_mf(mf_lambda);       //mf 1   Lambda
>     assem_boundary.push_mf(mf_U_line);       //mf 2   Ui
>     assem_boundary.push_mf(mf_interpolate);  //mf 3     LS
>     assem_boundary.push_mf(mf_U_interpolate);//mf 4     Ub
>     assem_boundary.push_data(LS_small);      //   data 1
>     assem_boundary.push_data(LAMBDA);        //   data 2
>     assem_boundary.push_vec(RHS_L_ui);       //   vec 1
>     assem_boundary.push_vec(RHS_L_ub);       //   vec 2
>     assem_boundary.push_vec(RHS_L_L);              //   vec 3
>     assem_boundary.assembly();
>
>
> I'm getting an error:
> ============================================
>
> |      An error has been detected !!!      |
>
> ============================================
> Error in getfem/bgeot_tensor.h, line 133 :
> index out of range
>
> The problem comes due to the assembly of V$2 (I checked that by 
> commenting out the other vectors). You can check that dimensions in 
> the assembly are correct (all the components of that assembly are 
> listed in the code). Is this a problem in the vectorized interpolated fem?
>
> Thank you again,
>                    Andriy
>
>
>
>
> -----Original Message-----
> From: Yves Renard [mailto:address@hidden
> Sent: Mon 1/26/2009 5:55 PM
> To: Andriy Andreykiv - 3ME
> Subject: Re: vectorial problem with interpolated fem
>
>
> Dear Andriy,
>
> The ser_qdim has to be done before the 
> mf_U_interpolate.set_finite_element().
>
> Yves.
>
> On Monday 26 January 2009 17:31, you wrote:
> > Dear Yves,
> >
> > I also tried to run the problem with getfem++3.1 (as opposed to 
> > 3.0.1) and I'm getting the same error message.
> >
> > Error in getfem_mesh_fem.cc, line 112 :
> > Incompatibility between Qdim=1 and target_dim 2 of FEM_UNKNOWN
> >
> > at line 3:
> >
> > 1:        getfem::mesh_fem mf_U_interpolate(*pline_mesh_deformed);
> > 2:        getfem::pfem
> > ifem_u=getfem::new_interpolated_fem(mf_bone_deformed,mim_line); 3:
> > mf_U_interpolate.set_finite_element(pline_mesh_deformed->convex_index(),i
> >fe m_u); 4:        mf_U_interpolate.set_qdim(pline_mesh_deformed->dim());
> >
> >
> > Best regards,
> >                  Andriy
> >
> >
> >
> >
> >
> >
> >
> > --------------------------------------------------------------------
> > ----
> > Andriy Andreykiv,
> > Postdoctoral researcher
> >
> > Leiden University Medical Center   |   Delft University of Technology
> > Biomechanics and Imaging Gr.       |   Fundamentals and Microsystems Gr.
> > Department of Orthopaedics         |   Faculty of 3mE
> > Room J09-127                       |   Room 8B-3-22
> > PO Box 9600                        |   Mekelweg 2
> > 2300 RC Leiden                     |   2628 CD Delft
> > The Netherlands                    |   The Netherlands
> > tel. +31-(0)71-5261566             |   tel. +31-(0)15-2786818
> > fax  +31-(0)71-5266743             |   fax  +31-(0)15-2789475
> > e-mail: address@hidden        |   e-mail: address@hidden
> > --------------------------------------------------------------------
> > ----
> >
> >
> >
> > -----Original Message-----
> > From: Yves Renard [mailto:address@hidden
> > Sent: Thu 1/22/2009 14:08
> > To: Andreykiv, A.
> > Subject: Re: vectorial problem with interpolated fem
> >
> >
> > Dear Andriy,
> >
> > A priori, the modification is not important to take into account 
> > vectorial fems (even if it can be optimized ...). I attached the 
> > modified file getfem_interpolated_fem.cc Can you test it and tell me 
> > if it works ? I am not completely sure of my modifications.
> >
> > Yves.
> >
> > On Tuesday 20 January 2009 14:10, you wrote:
> > > Dear Yves,
> > >
> > > A question on interpolated fem if I may. I need to use 
> > > interpolated fem with a vectorial problem, something that getfem 
> > > currently does not allow. Some time ago I was already asking you 
> > > about this and you suggested to reformulate my problem into a 
> > > scalar one, so I won't need the vector feature.
> > >
> > > I have a code now, which I'm, say, "enhancing" with interpolated 
> > > fem and reformulation of this code in scalar terms would be a 
> > > considerable time investment. Therefore, would you be kind to 
> > > briefly explain how difficult it is to develop a vectorial 
> > > interpolated fem and what steps need to be undertaken. Despite my 
> > > experience with getfem, I still have problems understanding the 
> > > internal structures, as the ones you have in getfem_interpolated_fem.h/cc.
> > >
> > > Is this, by any chance, anything you would be planning to include 
> > > in the next versions of getfem?
> > >
> > > Best regards,
> > >                 Andriy
> > >
> > > 
> > >-------------------------------------------------------------------
> > >----
> > >- Andriy Andreykiv,
> > > Postdoctoral researcher
> > >
> > > Leiden University Medical Center   |   Delft University of Technology
> > > Biomechanics and Imaging Gr.       |   Fundamentals and Microsystems
> > > Gr. Department of Orthopaedics         |   Faculty of 3mE
> > > Room J09-127                       |   Room 8B-3-22
> > > PO Box 9600                        |   Mekelweg 2
> > > 2300 RC Leiden                     |   2628 CD Delft
> > > The Netherlands                    |   The Netherlands
> > > tel. +31-(0)71-5261566             |   tel. +31-(0)15-2786818
> > > fax  +31-(0)71-5266743             |   fax  +31-(0)15-2789475
> > > e-mail: address@hidden        |   e-mail: address@hidden
> > > 
> > >-------------------------------------------------------------------
> > >----
> > >-
> > >
> > > Le mercredi 11 juillet 2007 19:14, Andriy Andreykiv a écrit :
> > > > Dear Yves,
> > > >
> > > >         What I want to do with interpolated fem is the 
> > > > following. I want to simulate a composite material (something 
> > > > like a car distribution belt), where the matrix of the composite 
> > > > is simulated with hexahedral elements while the fibers, that 
> > > > enforce the composite material are simulated with trusses. I 
> > > > want those trusses not to be conformal with the mesh of the 
> > > > matrix material, hence I want to impose a weak constraint that 
> > > > displacement of the fibers should be equal to the displacement 
> > > > of the matrix material in the location where those trusses 
> > > > fibers cross the matrix (this way I would attach the fibers 
> > > > inside the matrix  material). Hence, again this is conceptually 
> > > > similar to fictitious domain method for fluid-structure 
> > > > interection, where velocity of the fluid on the boundary of the 
> > > > structure should be equal the velocity of the structure and this 
> > > > constraint is imposed in a weak sence.
> > > >         I probably can reformulate my vectorial problem into a  
> > > > set of scalar problems and then use interpolated fem, and then 
> > > > assemble a vectorial problem again. Would that be a solution?
> > > >
> > > >I suppose you need the mass matrix beetween the two elements ?
> > > >  (because, if you only need a basic interpolation, in the sense 
> > > >of  Lagrange elements, you can have the interpolation matrix 
> > > >without  using interpolated fems).
> > > >
> > > >Yes, you can indeed try to reformulate your problem as a set of 
> > > >scalar problems and use the scalar interpolated fem. However, if 
> > > >you use only  scalar element (not intric vectorial ones), the 
> > > >adaptation of  interpolated fems to vectorial fems should not be 
> > > >to much  complicated. If you need, i can have a look to this when 
> > > >i will not  be too busy.
> > > >
> > > >
> > > >Yves.
> > > >
> > > >-----------------------------------------------------------------
> > > >-----
> > > >-- - Yves Renard (address@hidden)       tel : (33)
> > > > 04.72.43.80.11 Pole de Mathematiques,                       fax :
> > > > (33) 04.72.43.85.29 Institut Camille Jordan - CNRS UMR 5208
> > > >  INSA de Lyon, Universite de Lyon
> > > >  20, rue Albert Einstein
> > > >  69621 Villeurbanne Cedex, FRANCE
> > > >  http://math.univ-lyon1.fr/~renard
> > > >-----------------------------------------------------------------
> > > >-----
> > > >-- -

-- 

  Yves Renard (address@hidden)       tel : (33) 04.72.43.87.08
  Pole de Mathematiques, INSA-Lyon             fax : (33) 04.72.43.85.29
  20, rue Albert Einstein
  69621 Villeurbanne Cedex, FRANCE
  http://math.univ-lyon1.fr/~renard

---------