problem in coding FEM for 2D incompressible flow with triangle elemets
 

Dear all
I have written a FEM code for 2D incompressible flow, but I have encountered a problem. In my code every thing seems OK but the results are dreadful. I check everything (almost) and find nothing. I'm attaching all the functions I wrote. Please check them and let me know if there is problem or bug. It turns out to be my nightmare....:confused:

Bests
Mahdi

PS: just run main.m and see the results. You can check the results in TECPLOT but reading data file (fem.dat)

It would help if you described your application and method here.

Hi Jonas,

I used the method described in "nonlinear finite element analysis" by N.J. Reddy and it is a mixed method. I used linear triangle element for Ux and Uy and constant function for pressure (Np = 1). I'm writing this for lid driven problem to extend it for viscoelastic flow. The program is made up of five files, one main and four functions. read_data.m is reading the mesh information from node.txt, elements.txt and boundary conditions. KFlocal.m integrate local K and F using Gauss-quadrature rule. KFcalc.m calculate Kij and Fi for each integration point. I tried to write the routines as simple as possible so anyone could follow the code. If there is any question about the code... the problem is that when I try to solve the final equation set, [K]{U}={F}), MATLAB tells that [K] in badly scales or singular...

Do anybody have any idea about this problem?

I will not debug your code, and from the information you gave, it can be various things, and no one will be able to tell you exactly what is wrong.

Here's a few idea :

- Is you're mesh connectivity correct?
- Are every boundary conditions properly defined?
- Do you use a sufficient number of integration points?
- Are the element matrix (K) and vector (F) ok?
- Is the global assembly properly done?