[Getfem-users] Displacement fields in uniaxial tensile test

[Xin Liu]

Dear GetFEM-Users,

I wrote a simple 2D uniaxial tensile problem (see the mesh and coordinates in the attachment) mainly based on the example "thermo_elasticity_electrical_coupling".

The code is:

"

#include "getfem/getfem_mesh.h"

#include "getfem/getfem_mesh_fem.h"

#include "getfem/getfem_import.h"

#include "getfem/getfem_generic_assembly.h"

#include "gmm/gmm.h"

#include "getfem/getfem_model_solvers.h"

#include "getfem/getfem_mesher.h"

#include "getfem/getfem_fem.h"

using std::endl; using std::cout;

using bgeot::base_small_vector;

typedef getfem::model_real_plain_vector plain_vector;

int main(){

getfem::mesh m;

getfem::import_mesh("gmsh:./Simple.msh",m);

getfem::mesh_region border_faces;

getfem::outer_faces_of_mesh(m, border_faces);

getfem::mesh_region fb1 = getfem::select_faces_of_normal(m, border_faces,

base_small_vector( 1., 0.), 0.01);

getfem::mesh_region fb2 = getfem::select_faces_of_normal(m, border_faces,

base_small_vector( -1., 0.), 0.01);

m.region(1) = fb1; m.region(2) = fb2;

// Set finite element methods

getfem::pfem pf = getfem::fem_descriptor("FEM_QK(2,1)");

getfem::mesh_fem mfu(m); // Finite element for the elastic displacement

mfu.set_qdim(2);

mfu.set_finite_element(pf);

getfem::mesh_im mim(m); // Integration method

mim.set_integration_method(bgeot::dim_type(2));

getfem::model md;

md.add_fem_variable("u", mfu); // Displacement of the structure

// Membrane elastic deformation

double E = 21E6;

double nu = 0.3;

double F = 100E2;

double clambda = E*nu/((1+nu)*(1-2*nu)); // First Lame coefficient (N/cm^2)

double cmu = E/(2*(1+nu)); // Second Lame coefficient (N/cm^2)

double clambdastar = 2*clambda*cmu/(clambda+2*cmu); // Lame coefficient

md.add_initialized_scalar_data("cmu", cmu);

md.add_initialized_scalar_data("clambdastar", clambdastar);

getfem::add_isotropic_linearized_elasticity_brick

(md, mim, "u", "clambdastar", "cmu");

getfem::add_Dirichlet_condition_with_multipliers

(md, mim, "u", bgeot::dim_type(2-1), 2);

md.add_initialized_fixed_size_data("Fdata", base_small_vector(F,0.));

getfem::add_source_term_brick(md, mim, "u", "Fdata", 1);

gmm::iteration iter(1E-9, 1, 100);

getfem::standard_solve(md, iter);

//

// Solution export

plain_vector U(mfu.nb_dof());

gmm::copy(md.real_variable("u"), U);

bgeot::base_vector coeff;

getfem::slice_vector_on_basic_dof_of_element(mfu,U,1,coeff);

std::vector<std::string> v = { "Point a ", "Point b ", "Point c ", "Point d "};

cout << "Displacements for element 1: \n";

for (int i=0; i < coeff.size()/2; i++)

cout << v[i] << "u1: " << coeff[2*i] << "; u2： " << coeff[2*i+1] << endl;

return 0;

}

"

I output the displacement at points a, b, c, and d (see attachment).
Point a  u1: 0.000145356;  u2： 1.79962e-05
Point b  u1: 0.000161904;  u2： 6.42766e-05
Point c  u1: 0.000306614;  u2： 2.5169e-05
Point d  u1: 0.000313698;  u2： 7.23111e-05

I expect the u1 value of a and b is the same and c and d is the same in this problem since this is just an uniaxial tensile test. So does the difference come from numerical errors or some mistakes in the code?

Any comments/suggestions will be greatly appreicated.

Best regards,

--

======================================

Xin Liu

PhD Candidate, Graduate Research Assistant

Multiscale Structural Mechanics Group

https://sites.google.com/view/xinliu8809/

School of Aeronautics & Astronautics

Purdue University

West Lafayette, Indiana 47907-2045

======================================

SimpleStructure.png
Description: PNG image