AceFEM: Divergence in iterative procedure (Newton-Raphson) for fine meshes

Question

I am trying to model a problem of a nearly incompressible $10~\rm{m} \times 2~\rm{m}$ beam with a uniformly distributed end load. The beam has a Young's modulus of $200~\rm{Pa}$ and a Poisson's ratio of $\nu=0.49995$ and I am using a mixed displacement-pressure (Q1P0) formulation with a perturbed Lagrangian approach.

I am using meshes of $2^N \times 2^N$ Q1P0 elements and for $N\leq6$ to solution procedure is working fine. However, for $N>6$ the procedure diverges very early on. I have tried changing the initial load and the load step size in the solution procedure, and I have tried changing the material properties. Whatever I try I cannot get the solution procedure to work for $N>6$, which is a problem because I need solutions for at least up to $N=8$.
Note: The solution procedure works fine with standard Q2 elements for $N=8$ and I can't figure out why it doesn't work for the Q1P0 elements.

The solution procedure I am using is given below (unfortunately I don't know how to make the Greek characters display properly in the code block)

SMTAnalysis[];
tolNR = 10^-5; 
maxNR = 500; 
targetNR = 100;
\[Lambda]Max = 1; 
\[Lambda]0 = \[Lambda]Max/1000; 
\[CapitalDelta]\\[Lambda]Min = \[Lambda]Max/10000; 
\[CapitalDelta]\[Lambda]Max = \
\[Lambda]Max/100;

SMTNextStep["\[Lambda]" -> \[Lambda]0];
While[
  While[
   step = 
    SMTConvergence[tolNR, 
     maxNR, {"Adaptive BC", 
      targetNR, \[CapitalDelta]\[Lambda]Min, \
\[CapitalDelta]\[Lambda]Max, \[Lambda]Max}]
   , SMTNewtonIteration[];
   ];
  If[step[[4]] === "MinBound", SMTStatusReport["Analyze"]; 
   SMTStepBack[];];
  step[[3]] 
  , If[step[[1]], SMTStepBack[];];
  SMTNextStep["\[CapitalDelta]\[Lambda]" -> step[[2]]]
  ];

Plotting the displacement of the end of the beam against $N$ yields the following

It is clear that the behavior is good up until $N=6$.

The element code and notebook in which problem is solved (including solution procedure) are available here:
Q1P0 element
Convergence Test

Any help in improving the solution procedure would be appreciated!

Answer 1

I have noticed that your perturbation parameter goes to Infinity when Poisson's ratio goes to 0.5. Singularity always spells trouble for numerical procedures.

You have also formulated Lagrangian at the global level leading to the zeroes at the main diagonal. Additionally, close to singularity global matrix becomes ill-posed. The combined effects resulted in deteriorating convergence characteristics due to the MKL linear solver failure.

My advice is that you use static condensation of the Lagrangian (SMTNoDOFCondense option) to improve convergence.