[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: [Getfem-users] new_interpolated_fem problem
From: |
Andriy Andreykiv |
Subject: |
Re: [Getfem-users] new_interpolated_fem problem |
Date: |
Thu, 1 Mar 2007 17:48:43 +0100 |
User-agent: |
KMail/1.9.5 |
Dear Yves,
Thank you very much for your response. You were right and that helped.
I just
noticed that what I was trying to do was already done in
test_interpolated_fem.cc (in function test2) in the tests directory.
One more question if I may. After I received a mass matrix M which was
calculated from the shape fucntions of two overlapping meshes I have to
relate somehow the rows and columns of that matrix to the global degrees of
fredom of the original interpolated meshes.
If mq is the original mesh_fem (say some quad mesh) and m_intrpolate
is the
interpolation of mq on some mesh_fem, called ml (segment mesh that crosses
the quad mesh), then my mass matrix M has dimensions
M(m_interpolate.nb_dof(), ml.nb_dof). Imagine that I'm using this matrix as a
constraint matrix for an elliptical problem, where I'm imposing a Dirichlet
boundary condition in a weak sense as in a fictitious domain method, with
Lagrange multiplier:
| A B | { U }
| B^T 0 | {Lambda}
here U is a solution, defined on mq, A is an elliptic stiffness (mq.nb_dof X
mq.nb_dof) matrix, B is a constraint matrix and lambda is a Lagr.mult,
defined on the segment mesh_fem ml. Clearly B is bigger than M, as B has
dimentions B(mq.nb_dof, ml.nb_dof), so I have to relate rows of M to the
corresponding rows of B. I've looked at the code for the interpolated_fem,
but I didn't directly see how can I receive the index_of_dof for mq, that
correspond to the rows of M (which is probably numbered as the DOF's of
m_interpolate), so I can insert the rows of M into the correct places in B.
Would you be so kind to have a look at it?
Thank you very much in advance,
Andriy
On Wednesday 28 February 2007 15:38, you wrote:
> Le Mercredi 28 Février 2007 14:20, Andriy Andreykiv a écrit :
> > Good day gentlemen,
> >
> > I'm having problems with interpolating a finite element method from one
> > mesh to the other (getfem::new_interpolated_fem). What I want to do is
> > the following, below there is a code that first defines two quad elements
> > (in 2D space) and a line elemnet that crosses it arbitrarily in the
> > middle (there is a picture bellow). I want to calculate a mass matrix
> > (used for fictitious domain method) which is evaluated on that segment
> > element, but in fact is calculated as a product of shape functions from
> > the quads and the segment element (and the integration is over the
> > segment):
> >
> > M=\int_{segment} Base(the setment) X Base(the quads) dx
> >
> > I'm getting an error:
> >
> > interpolated_fem::update from context : convex_index = [0]
> > ============================================
> >
> > | An error has been detected !!! |
> >
> > ============================================
> > Error in ./gmm_opt.h, line 76 :
> > non invertible matrix
> >
> > I would really appreciate the help.
> > Best regards,
> > Andriy
>
> Hi Andriy,
>
> I see two errors in your program.
> The addition of the two quad element in the mesh should be writen:
>
> ind1[0]=quad_mesh.add_point(bgeot::base_node(0.0,0.0));
> ind1[1]=quad_mesh.add_point(bgeot::base_node(1.0,0.0));
> ind1[2]=quad_mesh.add_point(bgeot::base_node(0.0,1.0));
> ind1[3]=quad_mesh.add_point(bgeot::base_node(1.0,1.0));
> ind2[0]=ind1[1];
> ind2[1]=quad_mesh.add_point(bgeot::base_node(2.0,0.0));
> ind2[2]=ind1[3];
> ind2[3]=quad_mesh.add_point(bgeot::base_node(2.0,1.0));
> ind3[0]=line_mesh.add_point(bgeot::base_node(0.5,0.5));
> ind3[1]=line_mesh.add_point(bgeot::base_node(1.5,0.5));
> quad_mesh.add_parallelepiped(2,ind1.begin());
> quad_mesh.add_parallelepiped(2,ind2.begin());
>
> this particular numeration of the vertex come from the fact that
> quadrilaterons are defined as tensorial product of two segments in Getfem.
>
> Also, th line
> getfem::mesh_fem mq_interpolate(mq.linked_mesh());
>
> should be replaced by
> getfem::mesh_fem mq_interpolate(ml.linked_mesh());
>
> With this modifications, the program seems to work !
>
> Best regards,
>
> Yves.
>
> -------------------------------------------------------------------------
> Yves Renard (address@hidden) tel : (33) 04.72.43.80.11
> Pole de Mathematiques, INSA de Lyon fax : (33) 04.72.43.85.29
> Institut Camille Jordan - CNRS UMR 5208
> 20, rue Albert Einstein
> 69621 Villeurbanne Cedex, FRANCE
> http://math.univ-lyon1.fr/~renard
> -------------------------------------------------------------------------
--
__________________________________________
Andriy Andreykiv
Delft University of Technology
Faculty of Mechanical Engineering
Fundamentals of Microsystems
Mekelweg 2
2628 CD Delft
The Netherlands
E-mail : address@hidden
Tel. : (31) 15 2786818
Fax. : (31) 15 2789475
www : http://www-tm.wbmt.tudelft.nl/~andrico
private: (31) 6 47376804
- Re: [Getfem-users] new_interpolated_fem problem,
Andriy Andreykiv <=