Hi,
I am trying to use FreeFem++ to solve plasticity mechanics problems using an incremental iterative approach, with the specific code documentation being:
test.edp (15.6 KB)
tension.zip (307.4 KB) (it’s the mesh file)
, but I have discovered an issue:
When calculating stress, I use
Vh [u, v], [u1, v1];
Eh[int] sigma(3);
sigma[0] = E/(1-nu^2)*dx(u) + nu*E/(1-nu^2)*dy(v);
sigma[1] = nu*E/(1-nu^2)*dx(u) + E/(1-nu^2)*dy(v);
sigma[2] = 0.5*E/(1+nu)*(dx(v)+dy(u));
where u and v represent displacements, and sigma[3] represents stress.
The process of assigning values computed through dx(u) at integration points to the nodes of the sigma variable appears to introduce errors.
To verify this conjecture, I designed an extremely simple numerical example of linear elasticity
test_simple.edp (3.1 KB)
tension.zip (307.4 KB) (it’s the mesh file)
that should fully converge after the first iteration step. However, when I run the code:
solve A11([u,v],[u1,v1]) = int2d(Th)(epsilon(u,v)'*Copen(E,nu)*epsilon(u1,v1))
-int2d(Th)(epsilon(ubefore,vbefore)'*Copen(E,nu)*epsilon(u1,v1))
-int2d(Th)([sigma[0],sigma[1],sigma[2]]'*epsilon(u1,v1))
+on(1,u=0,v=-displacement)+on(2,v=displacement);
, it fails to converge. But when I run:
solve A11([u,v],[u1,v1]) = int2d(Th)(epsilon(u,v)'*Copen(E,nu)*epsilon(u1,v1))
-int2d(Th)(epsilon(ubefore,vbefore)'*Copen(E,nu)*epsilon(u1,v1))
-int2d(Th)(epsilon(ubefore,vbefore)'*Copen(E,nu)*epsilon(u1,v1))
+on(1,u=0,v=-displacement)+on(2,v=displacement);
it converges completely in the first step, just as expected.
It is strange that these two codes are essentially identical. I am unsure whether this discrepancy is caused by errors introduced when interpolating values from integration points to nodes. I would appreciate any help from friends here to resolve this issue.