[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
[Getfem-users] 2d plane strain/stress analysis
From: |
Umut Tabak |
Subject: |
[Getfem-users] 2d plane strain/stress analysis |
Date: |
Fri, 12 Mar 2010 14:03:17 +0100 |
User-agent: |
Mozilla-Thunderbird 2.0.0.22 (X11/20090707) |
Dear all,
I was wondering when using the
asm_stiffness_matrix_for_linear_elasticity function for a 2d problem(on
a mesh that is imported from gmsh), I would like to define a 2d
vectorial problem(a plane strain or stress one actually). So I set the
dimension to 2 since this is a 2s mesh. But when using this method in
the format given in the format of the manual
invalid data mesh fem (Qdim=1 required)
It has been sometime I have not been using the library, this can be a
simple mistake. The second question is that the thickness of the domain
in say z direction, is that taken as unit when computing the stiffness
matrices? And if I would like to define a different thickness, should I
use generic assembly procedure for the beam bending problems. I guess
there I can supply this as a parameter to the integration process.
Simple code is attached.
Best regards,
Umut
//
// 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_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 bgeot::scalar_type;
using bgeot::size_type;
//using namespace getfem;
//using namespace gmm;
using namespace std;
//
int main(int argc, char** argv)
{
try{
getfem::mesh mshGmsh;
import_mesh("gmshv2:./beamModel2d.msh", mshGmsh);
getfem::mesh_fem mf1(mshGmsh);
size_type dim = mshGmsh.dim();
std::cout << "Mesh dimension:" << dim << endl;
mf1.set_qdim(dim);
mf1.set_finite_element(mshGmsh.convex_index(),
getfem::fem_descriptor("FEM_QK(2, 1)"));
// check mesh dim and vectorial dimension equality
if(mshGmsh.dim() == mf1.get_qdim())
{
cout << "dimensions are equal" << endl;
cout << "number of dofs of the problem" << mf1.nb_dof() << endl;
}
// set integration methods
getfem::mesh_im intMethods(mshGmsh);
getfem::pintegration_method ppi =
getfem::int_method_descriptor("IM_QUAD(3)");
string im_name = getfem::name_of_int_method(ppi);
cout << im_name << endl;
intMethods.set_integration_method(ppi);
// define the matrix to use in the assembly process of
sparse_matrix SM(mf1.nb_dof(), mf1.nb_dof());
sparse_matrix MM(mf1.nb_dof(), mf1.nb_dof());
// define elastic constants
vector<double> lambda(mf1.nb_dof(), 70e9*0.33/((1+0.33)*(1-2*0.33)));
vector<double> mu(mf1.nb_dof(), 70e9/(2*(1+0.33)));
// assembly process
//asm_stiffness_matrix_for_linear_elasticity(SM, intMethods, mf1, mf1,
lambda, mu);
std::cout<<" --- stiffness matrix assembling --- " << std::endl;
//assem.assembly();
getfem::asm_stiffness_matrix_for_linear_elasticity(SM, intMethods, mf1,
mf1, lambda, mu);
std::cout<<" --- stiffness matrix assembled --- " << std::endl;
// std::cout<<" --- stiffness matrix: " << endl << SM << std::endl;
//std::cout<<" --- mass matrix assembling --- " << std::endl;
getfem::asm_mass_matrix(MM, intMethods, mf1, mf1);
std::cout<<" --- mass matrix assembled --- " << std::endl;
// write out the matrices in Harwell-Boeing format
string kFileName("Kmat.mm");
string mFileName("Mmat.mm");
// Harwell_Boeing_save(mFileName.c_str(), MM);
// Harwell_Boeing_save(kFileName.c_str(), SM);
// write out the matrices in Matrix-Market format
gmm::csc_matrix<double> SM2,MM2;
//gmm::copy(SM, SM2);
gmm::copy(MM, MM2);
MatrixMarket_save(mFileName.c_str(), MM2);
//MatrixMarket_save(kFileName.c_str(), SM2);
}
catch(std::exception &e){
std::cout << e.what() << std::endl;
}
return EXIT_SUCCESS;
}
- [Getfem-users] 2d plane strain/stress analysis,
Umut Tabak <=