# Power-law viscosity model for blood

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

 September 19, 2014, 05:21 Power-law viscosity model for blood #1 New Member   Filippo Join Date: May 2014 Posts: 25 Rep Power: 3 Hello friends, Hope everything is ok. I'm opening a new topic despite knowing that there are some other very similar. It happens that none clarifies me as I need... So, I'm using solver simpleFoam to simulate blood in a 2D geometry as a newtonian and non-newtonian fluid and the only model that is giving me problems is power -law. The power-law viscosity model has two parameters and, for blood, I found in the paper of Neofytou & Tsangaris (2006) the following values: k = 14.67×10−3 Pa.s^n n = 0.7755 I know that OpenFOAM uses kinematic viscosity and when I analysed the code of power-law model, noticed that the consistency index, k, as the same units of nu. So, I simply divided the value of k by rho (1040 kg/m3), and got a value of 1.411×10−5 m2/s for k. In a first phase, I used this value in my simulation and it diverged rapidly. In a second phase I have initialized all fields with potentialFoam and my simulation still doesn't converge. I noticed, in some tutorials, that k has values in the order of 2500 and I have some questions: Is anything escaping me? Am I supposed to obtain a value of k of this magnitude for blood with the data that I have available? One thing I do not understand is the n exponent in the units of k. Can somebody explain me why is this exponent there? What does it mean? I really need your help to proceed with my work. Sorry for my newbie question! Cheers, FS Last edited by fs.chemech; September 19, 2014 at 22:00.

 September 19, 2014, 14:40 #2 New Member   Filippo Join Date: May 2014 Posts: 25 Rep Power: 3 People, can someone give me a help on this question, please? It's really important to get my doubts solved! FS

 September 21, 2014, 13:21 #3 Super Moderator   Bruno Santos Join Date: Mar 2009 Location: Lisbon, Portugal Posts: 8,930 Blog Entries: 34 Rep Power: 92 Greetings FS, I quickly went to see your past posts and it seems I told you some months ago how important it is to increase the complexity in your test cases gradually, instead of jumping directly at the final problem Personally I never use the Power-law viscosity model, therefore I'm not familiar with it. If you can create and/or provide a simple example case, preferably with the expect analytical solution, for example something like the very first "cavity" tutorial explained in the User Guide, it should become pretty clear as to what is what . Best regards, Bruno __________________ OpenFOAM: Frequently Asked Questions | Useful links for building and using Forum: How to ask for help | Posting code and output with [CODE] My to-do list and when I'll be able to come to the forum: http://wyldckat.github.io And please: Read this before sending private messages to me

 September 23, 2014, 12:17 #4 New Member   Filippo Join Date: May 2014 Posts: 25 Rep Power: 3 Hi Bruno! First of all thank you for your answer! I remember well your advices! However, I think now the problem is slightly different and I think you could help me... I'm using a 2D T-shaped geometry (2 inlets and 1 outlet) and studying, with simpleFoam solver: 1) 3 levels of mesh refinement; 2) different schemes for divergence and pressure gradient; 3) Newtonian, Power-law and Bird-Carreau viscosity models. In a first phase, I studied the coarsest mesh: my simulations successfully converged when applying both type of schemes and the three viscosity models. For the intermediate and refined mesh, I only can attain convergence when applying newtonian viscosity model. I tried to manipulate relaxation factors for Power-law and Bird-Carreau but without success... Residuals for Ux and P stabilize far away from 1e-05! I'm using potentialFoam to initialize the fields in all my non-Newtonian simulations... Attached, there are several files for the intermediate refined mesh using a Power-law viscosity model, as well as a checkMesh file and my residuals evolution... Can you see them and give me some hints, please? Do you think that leaving simpleFoam and use nonNewtonianIcoFoam could be a solution? FS P.S.: as I said in my first message, I don't understand very well why there's an exponent "n" in the units of the consistency index "k"... In http://en.wikipedia.org/wiki/Power-law_fluid itit is possible to see that "K is the flow consistency index (SI units Pa•sn)"... Also, I don't know very well how this influences the value to input in powerLawCoeffs (transportProperties file)... checkMesh.txt p & U.txt transportProperties.txt controlDict fvScheme fvSolution.txt Residuals.gif

 September 24, 2014, 09:06 #5 Member   Florian Ries Join Date: Feb 2014 Location: Darmstadt, Germany Posts: 88 Rep Power: 3 Hi. the things about k, nu and n. When I have a look into /src/transportModels/incompressible/viscosityModels/powerLaw/powerLaw.C I see this: // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * // Foam::tmp Foam::viscosityModels:owerLaw::calcNu() const { return max ( nuMin_, min ( nuMax_, k_*pow ( max ( dimensionedScalar("one", dimTime, 1.0)*strainRate(), dimensionedScalar("VSMALL", dimless, VSMALL) ), n_.value() - scalar(1.0) ) ) ); } k_: flow consistency index n_: flow behaviour index strain-rate: strainRate() nu_: kin. viscosity k_ has the dimension of nu_. , because the strainRate() [1/sec] is multyplied by 1sec. ( "dimensionedScalar("one", dimTime, 1.0)"). The equation is like this: nu_ = k*{strain-Rate()*sec.}^(n-1) The flow behaviour index n is dimensionless. So everything is in the proper dimension. Is this the answer of your question?? kind regards Florian

 September 24, 2014, 09:21 #6 Member   Florian Ries Join Date: Feb 2014 Location: Darmstadt, Germany Posts: 88 Rep Power: 3 Hi, -----------fvSchemes--------- Why you use this: laplacianSchemes { default Gauss linear orthogonal; } snGradSchemes { default orthogonal; } instead of: laplacianSchemes { default Gauss linear corrected; } snGradSchemes { default corrected; } ???? Please post the complete case. Then we can see Kind regards Florian

 September 24, 2014, 10:12 #7 New Member   Filippo Join Date: May 2014 Posts: 25 Rep Power: 3 Florian, thank you very much for your answer!! Regarding the units of "k", I understood your explanation! Thank you very much! So, I'm using "orthogonal" and not "corrected" because my mesh is orthogonal. I saw (I don't remember where) that in these cases we can use this definition... Do you suggest me to change this? By the complete case, you mean my mesh and all the other folders ("0", "constant" and "system")? Cheers, FS

 September 24, 2014, 10:31 #8 Member   Florian Ries Join Date: Feb 2014 Location: Darmstadt, Germany Posts: 88 Rep Power: 3 Hi, Ah I see in your checkMesh-file: Mesh non-orthogonality Max: 0 average: 0 o.k. then everything is fine with that. yes with complete case I mean "0", "constant","system". Mesh is in PolyMesh. I don't need the blockMesh file. I will give it a shot at my PC. Then I will see what's going wrong with this case. kind regards Florian

 September 24, 2014, 11:34 #9 Member   Nadish Saini Join Date: Feb 2014 Location: Raleigh, North Carolina Posts: 37 Rep Power: 3 Hi 'fs.chemech', A comment on your original problem. I also faced a similar issue and could not figure out why the solution was diverging for non-newtonian law. Initalise your flow field with a newtonian model and then switch non-newtonian and it wont diverge. This solved the problem for me. Hope it helps!

 September 24, 2014, 12:18 #10 New Member   Filippo Join Date: May 2014 Posts: 25 Rep Power: 3 Florian, my compressed case (zip) is much bigger (~800 KB) than the max size allowed by the Forum. How can I send it to you? Nadish, thank you for your suggestion. You change to non-Newtonian model while Newtonian simulation is running?

 September 24, 2014, 12:26 #11 Member   Nadish Saini Join Date: Feb 2014 Location: Raleigh, North Carolina Posts: 37 Rep Power: 3 In the controlDict file you may set the option rumTimeModifiable to yes and then switch during runtime. But i guess the better option is to run the computation for say 0.1 sec, halt it, change the power law model and start from 0.1 sec. it will give a good initial value to the solution and i am pretty sure it wont diverge (given you have defined appropriate BC's)

 September 24, 2014, 13:02 #12 New Member   Filippo Join Date: May 2014 Posts: 25 Rep Power: 3 Nadish, I tried your suggestion... I know that with Newtonian model, my simulation converges in 290 iterations... So, I ran it until the 260ª iteration, stoped it, changed viscosity model to powerLaw and adapted startTime in controlDict file and resumed simulation again... Unfortunately, residuals Ux and p began to rise and stabilized again near 1e-04... Do you have any idea why this happens? Thank you for your time!

 September 24, 2014, 13:08 #13 Member   Nadish Saini Join Date: Feb 2014 Location: Raleigh, North Carolina Posts: 37 Rep Power: 3 Hi, You just need to initialize the values using Newtonian model. I don't think running upto 260 iterations is wise. Just run a few... Say 20... Iterations with Newtonian and then switch

 September 24, 2014, 13:13 #14 New Member   Filippo Join Date: May 2014 Posts: 25 Rep Power: 3 Nadish, I've done it... Residuals continue to stabilize around 1e-04...

 September 24, 2014, 13:36 #16 Member   Florian Ries Join Date: Feb 2014 Location: Darmstadt, Germany Posts: 88 Rep Power: 3 Hi again, if the fluctuations in your residuals will not decay, then you have transient effects in your flow. In this case your problem is unsteady and not steady. To get good avarage data, a steady simulation is not the best way. A better way is to track your field data over time and make an time-avarage (similar to LES or DNS). This can be very tricky, because the data are correlated over time (Autocorrelation). You need several integral time scale, until your results are good due to statistically accuracy. kind regards Florian

September 24, 2014, 13:40
#17
Member

Join Date: Feb 2014
Location: Raleigh, North Carolina
Posts: 37
Rep Power: 3
Quote:
 Originally Posted by fs.chemech Nadish, I've done it... Residuals continue to stabilize around 1e-04...

 September 24, 2014, 13:43 #18 New Member   Filippo Join Date: May 2014 Posts: 25 Rep Power: 3 Hi Florian, I don't have any experimental results to be based on... Yes, one of my suspicions is also that I have transient effects affecting my results... To solve my case as an unsteady case I have to change my solver to nonNewtonianIcoFoam isn't it? simpleFoam, doesn't run unsteady simulations, I think... You said that residuals are not good for checking accuracy! What, in your opinion, should I be based on? Cheers, FS

September 24, 2014, 13:52
#19
New Member

Filippo
Join Date: May 2014
Posts: 25
Rep Power: 3
Quote:
 Originally Posted by itchy "A better way is to track your field data over time and make an time-avarage (similar to LES or DNS). This can be very tricky, because the data are correlated over time (Autocorrelation). You need several integral time scale, until your results are good due to statistically accuracy."

Florian, can you explain me this a little better?
Sorry for my newbie questions...

 September 24, 2014, 14:22 #20 New Member   Filippo Join Date: May 2014 Posts: 25 Rep Power: 3 Hi Nadish, I wanted to say that I performed your suggestion. However, the simulation still does not converge ... Residuals continues to stabilize well above the desirable, around 1e-04! Thank you very much for your help!

 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 pbryant OpenFOAM Bugs 5 June 18, 2013 23:53 Qiaol618 Main CFD Forum 8 June 9, 2012 06:43 mannobot FLUENT 1 April 23, 2010 09:40 cpplabs OpenFOAM Running, Solving & CFD 1 February 13, 2008 09:09 Tim Phoenics 0 September 17, 2003 16:08

All times are GMT -4. The time now is 15:55.