Register Blogs Members List Search Today's Posts Mark Forums Read

July 17, 2012, 10:39
#1
Member

Elisabet Mas de les Valls
Join Date: Mar 2009
Location: Barcelona, Spain
Posts: 63
Rep Power: 9
Hi OpenFOAMers!
I've recently found a really strange behaviour on the calculation of the gradients on version 1.6-ext:

Starting from a uniform field, called P (partial pressure), which is uniform and equals 1, even at the boundaries, I want to evaluate the gradient of P, which is expected to be zero. However, it is not!!

Details on the case and code:

1.- the mesh is uniform, hence non-orthogonal corrections are not required. CheckMesh is o.k.
2.- boundary conditions of P are own made based on fixedValue type but, for this simple verification case, the value is set to 1 (and can be checked when looking at P values on the boundaries)
3.- two different strategies are used for the evaluation of the gradient:
3.a.
Code:
dPdx= fvc::reconstruct(sndPdx);
which I guess is better for non-orthogonal meshes
3.b.
Code:
In the attached figure you can find P map, dPdx_sn (from option 3.a.) and dPdx (from option 3.b.). As the mesh is orthogonal, both options yield to the same results. But the gradient is NOT zero!!!!

What's wrong? Am I going crazy?

any help is appreciated....

elisabet
Attached Images
 dPdx_small.jpg (35.6 KB, 14 views)

 July 17, 2012, 11:26 #2 Senior Member     Anton Kidess Join Date: May 2009 Location: Delft, Netherlands Posts: 1,141 Rep Power: 20 The gradient is (Phi_P - Phi_f)/mag(d). I'll assume d=1e-4m, and you get a small error in the computation of the difference of 1e-16 (that's about the precision you can get using double if I'm not mistaken). Dividing that by d leads to 1e-12, about the order of errors you are seeing. Further reading: http://en.wikipedia.org/wiki/Machine_epsilon __________________ *On twitter @akidTwit *Spend as much time formulating your questions as you expect people to spend on their answer. *Join the OpenFOAM stackexchange Q&A site: http://area51.stackexchange.com/prop...oHPxcPqde7HtA2 Last edited by akidess; July 17, 2012 at 11:28. Reason: Added machine epsilon link

 July 18, 2012, 07:38 #3 Member   Elisabet Mas de les Valls Join Date: Mar 2009 Location: Barcelona, Spain Posts: 63 Rep Power: 9 Thanks akidess! Obviously that was the point! sometimes, when one doesn't get the right solution begins to be suspicious of everything... sorry! now the bad new: I must have a mistake somewhere else in my code... let's see... elisabet

 July 18, 2012, 07:45 #4 Senior Member     Anton Kidess Join Date: May 2009 Location: Delft, Netherlands Posts: 1,141 Rep Power: 20 No worries and good hunting! __________________ *On twitter @akidTwit *Spend as much time formulating your questions as you expect people to spend on their answer. *Join the OpenFOAM stackexchange Q&A site: http://area51.stackexchange.com/prop...oHPxcPqde7HtA2

 Thread Tools Display Modes Linear Mode

 Posting Rules You may not post new threads You may not post replies You may not post attachments You may not edit your posts BB code is On Smilies are On [IMG] code is On HTML code is OffTrackbacks are On Pingbacks are On Refbacks are On Forum Rules

 Similar Threads Thread Thread Starter Forum Replies Last Post Chander CFX 5 January 7, 2015 04:00 senD Fluent UDF and Scheme Programming 9 September 18, 2014 07:29 fruitkiwi Main CFD Forum 4 June 12, 2012 01:12 rks171 Main CFD Forum 1 January 4, 2012 17:54 Raman Chadha FLUENT 0 January 12, 2005 19:00

All times are GMT -4. The time now is 09:56.