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
--
//++++++++++++++++++++++THE CODE++++++++++++++++++++++++
#include <getfem_model_solvers.h>
#include <gmm.h>
#include <getfem_interpolated_fem.h>
/* some Getfem++ types that we will be using */
using bgeot::base_small_vector; /* special class for small (dim<16) vectors */
using bgeot::base_node; /* geometrical nodes(derived from base_small_vector)*/
using bgeot::scalar_type; /* = double */
using bgeot::size_type; /* = unsigned long */
using bgeot::base_matrix; /* small dense matrix. */
/* definition of some matrix/vector types.
* default types of getfem_model_solvers.h
*/
typedef getfem::modeling_standard_sparse_vector sparse_vector;
typedef getfem::modeling_standard_sparse_matrix sparse_matrix;
typedef getfem::modeling_standard_plain_vector plain_vector;
int main(int argc, char *argv[]) {
try {
// Building a mesh of two quads and a segment that crosses them in the middle
// ___________________________________
// | | |
// | | |
// | | |
// | +-----------------------------+ |
// | | |
// | | |
// | | |
// |_________________|__________________|
//
getfem::mesh quad_mesh;
getfem::mesh line_mesh;
std::vector<bgeot::size_type> ind1(4);
std::vector<bgeot::size_type> ind2(4);
std::vector<bgeot::size_type> ind3(2);
ind1[0]=quad_mesh.add_point(bgeot::base_node(0.0,0.0));
ind1[1]=quad_mesh.add_point(bgeot::base_node(0.0,1.0));
ind1[2]=quad_mesh.add_point(bgeot::base_node(1.0,1.0));
ind1[3]=quad_mesh.add_point(bgeot::base_node(1.0,0.0));
ind2[0]=ind1[2];
ind2[1]=ind1[3];
ind2[2]=quad_mesh.add_point(bgeot::base_node(2.0,0.0));
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());
line_mesh.add_convex(bgeot::simplex_geotrans(1,1),ind3.begin());
// mesh_fem for the two quad elements
getfem::mesh_fem mq(quad_mesh);
mq.set_finite_element(mq.linked_mesh().convex_index(),getfem::fem_descriptor("FEM_QK(2,1)"));
mq.set_qdim(1);
// mesh_fem and integration method for the segment element
getfem::mesh_fem ml(line_mesh);
ml.set_finite_element(getfem::fem_descriptor("FEM_PK(1,1)"));
ml.set_qdim(1);
getfem::mesh_im mim(line_mesh);
getfem::pintegration_method
ppi=getfem::int_method_descriptor("IM_GAUSS1D(3)");
mim.set_integration_method(ml.linked_mesh().convex_index(),ppi);
// CREATION OF THE INTERPOLATION MEHS_FEM, THAT IS SUPPOSED TO INTERPOLATE
// THE QUAD'S SHAPE FUNCTIONS ON THE SEGMENT
getfem::mesh_fem mq_interpolate(mq.linked_mesh());
getfem::pfem ifem=getfem::new_interpolated_fem(mq,mim); //THIS DOESN'T
// WORK, DON'T KNOW WHY !!!!!!
mq_interpolate.set_finite_element(ifem);
// mass matrix on the segment created from the shape function of the quads
and the segment
sparse_matrix M(6,2);
getfem::asm_mass_matrix(M,mim,mq_interpolate,ml);
//
std::cout<<" M="<<M<<std::endl;
//
}
DAL_STANDARD_CATCH_ERROR;
return 0;
}
_______________________________________________
Getfem-users mailing list
[email protected]
https://mail.gna.org/listinfo/getfem-users
