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

msalamih · May 24, 2014, 13:48

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....

Bests
Mahdi

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

Jonas Holdeman · May 25, 2014, 11:37

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

msalamih · May 25, 2014, 13:24

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...

msalamih · May 26, 2014, 16:22

Do anybody have any idea about this problem?

Alex C. · May 26, 2014, 20:19

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?

Thread Tools	Search this Thread
Show Printable Version Email this Page	Search this Thread: Advanced Search
Display Modes
Linear Mode Switch to Hybrid Mode Switch to Threaded Mode

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
Complex problem: Open Channel Flow	lett	FLUENT	18	November 10, 2018 18:30
transient, impregnating flow problem	fgommer	FLUENT	0	February 29, 2012 16:10
High Reynolds Number Incompressible Flow	Nick Lavery	Main CFD Forum	3	February 19, 2002 15:47
fluid flow fundas	ram	Main CFD Forum	5	June 17, 2000 21:31
How good is CFD?	kai	Main CFD Forum	39	April 7, 2000 12:48

May 25, 2014, 11:37		#2
Jonas Holdeman Senior Member Jonas T. Holdeman, Jr. Join Date: Mar 2009 Location: Knoxville, Tennessee Posts: 128 Rep Power: 18	It would help if you described your application and method here.

May 25, 2014, 13:24		#3
msalamih New Member Mahdi Salami Join Date: Feb 2014 Posts: 6 Rep Power: 12	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...

May 26, 2014, 16:22		#4
msalamih New Member Mahdi Salami Join Date: Feb 2014 Posts: 6 Rep Power: 12	Do anybody have any idea about this problem?

May 26, 2014, 20:19		#5
Alex C. Member Join Date: Jul 2013 Posts: 55 Rep Power: 12	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?