[Top][All Lists]

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: [Getfem-users] different mesh_fem definitions for different regions

From: Umut Tabak
Subject: Re: [Getfem-users] different mesh_fem definitions for different regions
Date: Thu, 30 Jul 2009 14:09:18 +0200
User-agent: Thunderbird (X11/20090409)

Renard Yves wrote:

Umut Tabak <address@hidden> a écrit :

Renard Yves wrote:

Dear Umut,

If your region is a part of a mesh (not a boundary)
Hi again,

I tried to apply a sample analysis where I can assemble the fluid part(scalar) and the structural part(displacement). Now I have a problem concerning the coupling matrix assembly. I attached the code and the necessary msh files. Now the questions are(better to answer these after taking a look at the source file attached. In brief, I have 3 mesh structures and 3 mesh_fems attached to these mesh structures and the operations should be clear for you.)

+ To have this kind of analysis(coupled structural-acoustic), we have to set the mesh_fem dimensions to 3 and 1 for the structural and fluid parts respectively, if there is only one mesh_fem definition, how can I assign different dimensions to different regions of my mesh?

+ The second question follows accordingly, how to assign different integration methods to different regions of the mesh? Assembly can be done on regions as you suggested before.

+ Without the definition of separate mesh_fem objects, how can I allocate memory for the matrices that I am going to use for the assembly process? I mean to get the specific dof numbers associated to regions.

+ Lastly, with the attached file, I can assemble structural domain and apply boundary conditions on that, also the fluid part, but i have a problem related to a term like

 \bold{K}_{up}=\int_{S}{\bold{N}_{s}}^{T} \bold{n}{\bold{N}_{p}}dS

that I would like to assemble on the interface of the two above defined domains.

Shape functions N_s will come from the structural domain and N_p will be the contributions of the fluid domain. But i could not assemble this term over the surface defined in the mesh, I understand the error message but I could not think a way to fix it. The error is ''error: the mesh_fem/mesh_im live on different meshes!''.

I would appreciate further help on this.

Thanks and best regards,

// standard headers
#include <iostream>
#include <vector>
#include <string>
#include <cstdlib>
#include <cassert>
#include <stdexcept>
// gmm++ headers, pay attention that this is a template
// library, so that no build is required.
#include <gmm/gmm.h>
#include <gmm/gmm_inoutput.h>
// getFem++ headers 
#include <getfem/getfem_import.h>
#include <getfem/dal_bit_vector.h>
#include <getfem/getfem_mesh.h>
#include <getfem/getfem_mesh_fem.h>
#include <getfem/getfem_partial_mesh_fem.h>
#include <getfem/getfem_mesh_im.h>
#include <getfem/getfem_integration.h>
#include <getfem/getfem_assembling.h>
#include <getfem/getfem_modeling.h>
#include <getfem/getfem_regular_meshes.h>
typedef getfem::modeling_standard_sparse_matrix sparse_matrix;
//using namespace getfem;
//using namespace gmm;
using namespace std;
using bgeot::scalar_type;
int main(int argc, char** argv)
    // if(argc!=3){
    //   throw std::runtime_error("Error in input parameters\nUsage: <binary> 
<filename.msh> <gmshFormatSpecifier>");
    // }
    getfem::mesh mshGmsh1;
    getfem::mesh mshGmsh2;
    getfem::mesh mshGmsh3;
    //import_mesh(argv[1], argv[2], mshGmsh);
    getfem::import_mesh("gmshv2:./boxSolidDom1.msh", mshGmsh1);
    getfem::import_mesh("gmshv2:./boxSolidDom2.msh", mshGmsh2);
    getfem::import_mesh("gmshv2:./boxSolidDom3.msh", mshGmsh3);
                // optional print out for the convexes(elements) of a specific
                // region
                cout << "Number of elements in region 1 is " << 
mshGmsh1.nb_convex() << endl;
                cout << "Number of elements in region 2 is " << 
mshGmsh2.nb_convex() << endl;
                cout << "Number of elements in region 3 is " << 
mshGmsh3.nb_convex() << endl;
    // can also be an array of vectors
                // dal::bit_vector elemsInRegion1;
                // dal::bit_vector elemsInRegion2;
                // dal::bit_vector elemsInRegion3;
                // element definitions 
                // for (getfem::mr_visitor i(region1); !i.finished(); ++i) 
                // for (getfem::mr_visitor i(region2); !i.finished(); ++i) 
                // for (getfem::mr_visitor i(region3); !i.finished(); ++i) 
    // assign fem definitions to the regions
                getfem::mesh_fem meshFemDom1(mshGmsh1);
                getfem::mesh_fem meshFemDom2(mshGmsh2);
                getfem::mesh_fem meshFemDom3(mshGmsh3);
                // define also a mesh fem for the right hand side 
                // for boundary condition application
                getfem::mesh_fem rhsAssembly(mshGmsh1);
                getfem::mesh_fem lameConsts(mshGmsh1);
                // structural domain 
meshFemDom1.set_finite_element(getfem::fem_descriptor("FEM_QK(3, 1)"));
    cout << "Number of dofs of the structural domain: " << meshFemDom1.nb_dof() 
<< endl;
                                   getfem::fem_descriptor("FEM_QK(3, 1)"));
                           getfem::fem_descriptor("FEM_QK(3, 0)"));
                // set integration methods of different domains
    getfem::mesh_im intMethods1(mshGmsh1);
    getfem::pintegration_method ppi1 = 
                // define the matrix to use in the assembly process of
    sparse_matrix SM(meshFemDom1.nb_dof(), meshFemDom1.nb_dof());
    sparse_matrix MM(meshFemDom1.nb_dof(), meshFemDom1.nb_dof());
    // define the source term
    vector<scalar_type> B(meshFemDom1.nb_dof());
    // define elastic constants
    vector<double> lambda(lameConsts.nb_dof(), 70e9*0.33/((1+0.33)*(1-2*0.33)));
    vector<double> mu(lameConsts.nb_dof(), 70e9/(2*(1+0.33)));
    // bc values array
    vector<scalar_type> V(meshFemDom1.nb_dof(),0.0);
    V[7] = 10;
    cout << "number of dofs of the RHS term" << rhsAssembly.nb_dof() << endl;
    // assembly process
    std::cout<<" ---  stiffness matrix assembling --- " << std::endl;
    //getfem::asm_stiffness_matrix_for_linear_elasticity(SM, intMethods, 
meshFemDom1, rhsAssembly, lambda, mu);
    getfem::asm_stiffness_matrix_for_linear_elasticity(SM, intMethods1, 
meshFemDom1, lameConsts, lambda, mu);
    std::cout<<" ---  stiffness matrix assembled  --- " << std::endl;
    std::cout<<" ---       mass matrix assembling --- " << std::endl;
    getfem::asm_mass_matrix(MM, intMethods1, meshFemDom1,meshFemDom1);
    std::cout<<" ---       mass matrix assembled  --- " << std::endl;
    std::cout<<" --- source term assembling --- " << std::endl;
    getfem::asm_source_term(B, intMethods1, meshFemDom1, rhsAssembly, V, 2);
    //getfem::asm_source_term(B, intMethods, meshFemDom1, lameConsts, V, 2);
    std::cout<<" ---  source term assembled --- " << std::endl;
    // //
    std::cout<<" --- applying the boundary conditions --- " << std::endl;
    getfem::assembling_Dirichlet_condition(SM, B, meshFemDom1, 2, V);
    getfem::assembling_Dirichlet_condition(MM, B, meshFemDom1, 2, V);
    std::cout<<" --- applied the boundary conditions  --- " << std::endl;
    // dump matrices onto files
                gmm::csc_matrix<double> SM2,MM2;
    string kFileName("");
    string mFileName("");
    gmm::copy(SM, SM2);
    gmm::copy(MM, MM2);
    MatrixMarket_save(mFileName.c_str(), MM2);
    MatrixMarket_save(kFileName.c_str(), SM2);
    // fluid domain, Domain 2, scalar fem, only pressures are considered
                                   getfem::fem_descriptor("FEM_QK(3, 1)"));
                // set integration methods of different domains
    getfem::mesh_im intMethods2(mshGmsh2);
    getfem::pintegration_method ppi2 = 
                sparse_matrix SMF(meshFemDom2.nb_dof(), meshFemDom2.nb_dof());
    sparse_matrix MMF(meshFemDom2.nb_dof(), meshFemDom2.nb_dof());
          // assemble fluid only part
    getfem::generic_assembly assem;
    std::cout<<" ---  Fluid stiffness matrix assembling --- " << std::endl;
    std::cout<<" ---  Fluid stiffness matrix assembled  --- " << std::endl;
    // std::cout<<" ---  stiffness matrix:  " << endl << SM << std::endl;
    std::cout<<" ---  Fluid mass matrix assembling --- " << std::endl;
    getfem::asm_mass_matrix(MMF, intMethods2, meshFemDom2, meshFemDom2);
    std::cout<<" ---  Fluid mass matrix assembled  --- " << std::endl;
    // dump matrices of the fluid part
                gmm::csc_matrix<double> SM3,MM3;
    gmm::copy(SMF, SM3);
    gmm::copy(MMF, MM3);
    string kFileNameFluid("");
    string mFileNameFluid("");
    MatrixMarket_save(mFileNameFluid.c_str(), MM3);
    MatrixMarket_save(kFileNameFluid.c_str(), SM3);
                // set integration methods for the coupling surface 
    getfem::mesh_im intMethods3(mshGmsh3);
    getfem::pintegration_method ppi3 = 
    // allocate the coupling matrix
                sparse_matrix SMCoupling(gmm::mat_nrows(SM), 
    //  assemble coupling terms
    getfem::generic_assembly assemCoupling;
                assemCoupling.push_mf(meshFemDom1); // structural domain #1
    assemCoupling.push_mf(meshFemDom2); // fluid domain      #2
    std::cout<<" ---  Coupling matrix assembling --- " << std::endl;
    std::cout<<" ---  Coupling matrix assembled  --- " << std::endl;
    // dump matrices of the fluid part
                gmm::csc_matrix<double> SMC;
    gmm::copy(SMCoupling, SMC);
    string kFileNameCoupling("");
  catch(std::exception &e){
    std::cout << e.what() << std::endl;
  return EXIT_SUCCESS;


Attachment: boxSolidDom1.msh
Description: Mesh model

Attachment: boxSolidDom2.msh
Description: Mesh model

Attachment: boxSolidDom3.msh
Description: Mesh model

reply via email to

[Prev in Thread] Current Thread [Next in Thread]