A strange problem for a nonlinear fea

[ricewjx]

I am developing a nonlinear finite element program for the response of ferroelectric materials. Due to the electro-mechanical coupling, the nonlinearity is material nonlinear. I derived the consistent tangent stiffness matrix and are trying to test it. Three tests were carried out on a rectangular plate with both linear stiffness matrix and consistent tangent matrix and the results are as following:

1. Uniform tension without polarization (electrical term). The results are what we expected for both kinds of stiffness matrics: homogenous.

2. Still uniform tension but with the material initially poled. The results are also uniform although the linear stiffness matix will take more steps to converge. This is right since we expect the consistent tangent stiffness will increase the convergent speed.

3. Uniform electrical loading. Load control was used in 1st and 2nd test. Due to the formulation of the finite element method, the 3rd test is using displacement control. Here, the displacements are fixed to zero and the electrical degree of freedom is increased step by step. Now the funny thing happens. For the linear stiffness case, the results are homogeneous even with very small or large step size. For the consistent stiffness case, it will converge but the results are not homogeneous if a lot of steps are taken.

I was think of what might cause the problem. At first I thought it maybe my derivation of the tangent stiffness has some mistakes. However, if that is the case, I would not expect the 2nd test to have correct solution. One possibility is the consistent tangent matrix is slightly unsymmetric but our solver is for symmetric matrix. So I made the matrix symmetric by (A+AT)/2. I don't know whether this will cause the problem. Another possibility is the numerical error but I don't know how I can check it.

Can anyone tell me what might cause the problem and how I can do to eliminate the problem? Thanks a lot for your help!

 

[EricZhao]

You'd better use some fumulars to explain what you did. Some terms you used might be spedific to your project. Few people could understand.

Since you talked about the treatment: (A+AT)/2. If your matrix:A has very small eigenvalue, a change of any entry of the matrix could result great error of the matrix inverse. If your matrix is supposed to be unsymmetric, don't make it symmetric anyway.

 

[ricewjx]

Thank you EricZhao! I will try to demonstrate the problem. I hope I can make it clear.
The consitutive relationships are plasticity like:

?=C?-?r)-hD-Pr) (1)
E=-h?-?r)+?D-Pr) (2)

?-stress, E-electric field, ?-strain, ?r-remanent strain, D-polarization, Pr-remanent polarization, C, h, ? are material coefficients and are dependent on Pr.

The finite element formulation is:

Integrate(?:??+E:?D)dV=Integrate(t:?u-??(D:n))dS (3)

t-traction, ?-electric potential, n-unit vector normal to the surface.

The right handside of (3) is the external force F, and the left handside of (3) is the internal force R. For every increment of displacement, there's an error of ?R=R-F. So we have

K:?U=?R

Here K is the consistent tangent stiffness which involves with a lot of derivations, U is the nodal degree of freedom which includes displacement and electrical function ?. D=-partial ?/partial x. x is the coordinates

So for the 1st test, there is no Pr or electrical terms invovled and it's just like traditional plasticity: with the increase of traction, the material response will change from elasticity to plasticity. In the 2nd test, the material initially has a remanent polarization Pr0. The response will be nonlinear from the first step and electrical terms will play a role in it. For this two tests, the results are what we expected: homogeneous and consistent matrix will need less steps to converge than linear stiffness matrix. Here linear stiffness matrix means C, h, and ?.

In the 3rd test, instead increasing the traction, I increase ? to get uniform fields. With linear stiffness matrix, the result is homogeneous no matter how small the step size is. For the consistent stiffness matrix, the results will converge to inhomogeneous solution though it will have correct solution for the first several nonlinear response.

As I artificially made the K symmetric, this maybe the reason as EricZhao mentioned. Our finite element storage is skyline scheme. So if anyone can tell me where I can find a routin that solves nonsymmetric matrix with skyline storage, I will be really grateful.

 

[ricewjx]

Sorry. The greek characters can not be displayed. I have to make changes and post it again.
The consitutive relationships are plasticity like:

sigma=Ceps-epsr)-hD-Pr) (1)
E=-heps-epsr)+bD-Pr) (2)

sigma-stress, E-electric field, eps-strain, epsr-remanent strain, D-polarization, Pr-remanent polarization, C, h, b are material coefficients and are dependent on Pr.

The finite element formulation is:

Integrate(sigma:delta(eps)+E:delta(D))dV=Integrate(t:delta(u)-phi delta(D:n))dS (3)

t-traction, phi-electric potential, n-unit vector normal to the surface.

The right handside of (3) is the external force F, and the left handside of (3) is the internal force R. For every increment of displacement, there's an error of dR=R-F. So we have

K:dU=dR

Here K is the consistent tangent stiffness which involves with a lot of derivations, U is the nodal degree of freedom which includes displacement and electrical function psi. D=-partial psi/partial x. x is the coordinates

So for the 1st test, there is no Pr or electrical terms invovled and it's just like traditional plasticity: with the increase of traction, the material response will change from elasticity to plasticity. In the 2nd test, the material initially has a remanent polarization Pr0. The response will be nonlinear from the first step and electrical terms will play a role in it. For this two tests, the results are what we expected: homogeneous and consistent matrix will need less steps to converge than linear stiffness matrix. Here linear stiffness matrix means C, h, and b.

In the 3rd test, instead increasing the traction, I increase psi to get uniform fields. With linear stiffness matrix, the result is homogeneous no matter how small the step size is. For the consistent stiffness matrix, the results will converge to inhomogeneous solution though it will have correct solution for the first several nonlinear response.

As I artificially made the K symmetric, this maybe the reason as EricZhao mentioned. Our finite element storage is skyline scheme. So if anyone can tell me where I can find a routin that solves nonsymmetric matrix with skyline storage, I will be really grateful.