[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: |
Mon, 03 Aug 2009 14:52:24 +0200 |
User-agent: |
Thunderbird 2.0.0.21 (X11/20090409) |
Renard Yves wrote:
Umut Tabak <address@hidden> a écrit :
Renard Yves wrote:
Umut Tabak <address@hidden> a écrit :
Renard Yves wrote:
Dear Umut,
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?
You cannot ! Even if you have only one mesh representing both the
structure and the fluid, the best is to have two mesh_fem, each
defined only on the corresponding region. This is the easiest for
computing the coupling terms but of course it means that the meshes of
the structure and the fluid match to each other.
If you have two separate meshes you can try to compute the coupling
term using the object "interpolated fem" which interpolate a fem from
a mesh to another (see for instance tests/test_interpolated_fem.cc)
Dear Professor Renard,
That is the question I am still puzzling with, how to assign two
mesh_fems to different regions of a mesh, I am sorry again but let me
ask one more, I can iterate over regions as you have suggested, to
assign pfem definitions, but how to change the dimensions of the these
regions? As far as, I can understand from your definitions, I can have
one mesh and different mesh_fem definitions for different regions. I
could partly accomplish this task by assigning the pfem by iteration but
I still could not understand how to change the target dimension, Or do I
have a conceptual misunderstanding on this point? Or there is a trick to
do this that I am still missing? There is a simple code below to explain
my purposes.
Best regards,
Umut
A very simple code vhere I think the dimensions also be changed, but you
are the creator of the library.
// 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;
//
int main(int argc, char** argv)
{
try{
getfem::mesh mshGmsh;
getfem::import_mesh("gmshv2:./boxSolidCheck.msh", mshGmsh);
// region 1 is the fluid part of the model mesh definition
if(mshGmsh.has_region(1))
cout << "Region 1 is defined" << endl;
if(mshGmsh.has_region(2))
cout << "Region 2 is defined" << endl;
if(mshGmsh.has_region(3))
cout << "Region 3 is defined" << endl;
getfem::mesh_region region1(mshGmsh.region(1));
getfem::mesh_region region2(mshGmsh.region(2));
getfem::mesh_region region3(mshGmsh.region(3));
// optional print out for the convexes(elements) of a specific
// region
cout << "Number of elements in region 1 is " <<
region1.nb_convex() << endl;
cout << "Number of elements in region 2 is " <<
region2.nb_convex() << endl;
cout << "Number of elements in region 3 is " <<
region3.nb_convex() << endl;
// can also be an array of vectors
if(region1.is_only_convexes())
cout << "region 1 only contains convexes" << endl;
else
cout << "region 1 is not composed of convexes only" << endl;
//
dal::bit_vector elemsInRegion1= region1.index();
dal::bit_vector elemsInRegion2= region2.index();
dal::bit_vector elemsInRegion3= region3.index();
//
for(getfem::mr_visitor i(region1);!i.finished();++i){
cout << "Element :" << i.cv() << " in region1." << endl;
}
// assign fem definitions to the regions
getfem::mesh_fem meshFemCore(mshGmsh);
// region 1, scalar domain, pressure field, dimension 1
for(getfem::mr_visitor i(region1); !i.finished(); ++i){
meshFemCore.set_finite_element(i.cv(),getfem::fem_descriptor("FEM_QK(3,
1)"));
// should not dimensions set here??
}
// region 2, vectorial domain, displacement field, dimension 3
for(getfem::mr_visitor i(region2); !i.finished(); ++i){
meshFemCore.set_finite_element(i.cv(),getfem::fem_descriptor("FEM_QK(3,
1)"));
// should not dimensions set here??
}
}
catch(std::exception &e){
std::cout << e.what() << std::endl;
}
return EXIT_SUCCESS;
}