Finite Element Code: Convergence issues

Question

I am trying to code a simple FEM problem. A bar is fixed at the bottom edge, and a displacement is applied at the top. I want to apply the displacement gradually , so I have divided it into a number of steps. A rough outline of the code looks like this:

 %********************************************************* %          
> Newton-Raphson Control parameters
> %*********************************************************
> 
> uload=5;         % the total displacement  
  numu=10;         % number of steps 
  du=uload/numu;   % step size 
  tu=0;            %initial value
> tol=1e-2; 


>   for ii=1:numu  %Loop over total number of steps
> 
>       iter=0;
>       tu=tu+du;  %increment displacement at each step
> 

>          while(error>=tol)  %Iteration
>            iter=iter+1
>            Ku=sparse(sdof,sdof);    %Initialize stiffness matrix
>            [Ku] = Kmat(Ku);  %Create stiffness matrix by user defined function
>               Ru = sparse(sdof,1);  %Initialize residual force vector
>            [Ru] = Rvec(Ru,u);  %Calculate residual force vector by user defined function(consists of internal force only)
>            [Ku,Ru] = applyBC(Ku,Ru); %Apply BCs by user defined function

>            delu = sparse(sdof,1) ;  %Initialize solution vector       
>            delu = Ku\Ru ;   % Calculate displacement
             u=u+delu; %   %Increase displacement value
>            error=norm(Ru)  %Calculate Error (This step is probably causing the divergence due to selection of improper criterion)



>          end %iteraion  
> end %step

This code won't compile in MATLAB because all the functions are user defined and the input data is not provided. I just wanted to portray a rough outline of what I am trying to do.

This algorithm, works fine if I only have a few elements(suppose 100) defined. If I decrease the mesh size, the solutions tends to diverge, the error keeps on increasing. I figure that is probably because the way I calculate the error is not useful for this simple case.

Now my question is what can be the best convergence criterion for this problem? Since the problem is simple linear elastic, no non-linearity involved, I think the values calculated in the first iteration is correct, so no need for further iterations. But later I would like to extend this code to perform non-linear static analysis, so is there a convergence criterion that can be used for both? If not, then what can be used for each case?

P.S: I have modeled a similar problem in a commercial code, where obviously it works, but when I check the details of the solution, it says that it has used one iteration per step(which is expected). In my code with more elements, the iterations keep increasing after the 2nd step. Thus my doubt regarding the choice of convergence criterion.

Answer 1

I think your convergence tolerance is not enough. Since you obtain first few solutions with low convergence requirements you are moving away from the real results for the next solution steps. From the code, it seems you are using unrealistic values. In example, total displacement is 5. The tolerance is 0.01, they are crude numbers. So it is not easy to determine which values must be used. Without proper use of realistic values, it is not easy for the model to converge. For a nonlinear problem, suppose you find disp1=0.01 with corresponding force1=0.01, then for disp2=0.02, assume an acceptable solution for force2=0.019, but if without properly adjusting tolerance, you may find force2=0.018, then the force3 result will definitely change based on the force2 value you obtained.