[eberberovic]

Dear all,

I would like to share the implementation of a new boundary condition of dynamic contact angle for interFoam aplication. It is based on the investigations of Hoffman, firstly introduced by Kistler(1993) through a Hoffman function. Details can be found in:

@InCollection{Kistler93,
author = {Stephan F. Kistler},
title = {Hydrodynymics of Wetting},
booktitle = {Wettability},
pages = {311--429},
publisher = {Marcel Dekker Inc.},
year = {1993},
editor = {John C. Berg},
chapter = {6},
}

Basically, one of the problems in assigning an expression for the dynamic contact angle (in addition to modelling the contact line velocity) seems to be the appearance of a static contact angle. The static contact angle shows a hysteresis, and is therefore not single valued. In order to overcome this difficulty, a term is introduced in the function for the dynamic contact angle, through which the influence of the static contact angle is absorbed. The function is of the form, Kistler(1993):

thetaD = F_hoff(Ca + F_hoff^-1(theta0)),

where the symbols are :
thetaD – dynamic contact angle,
F_hoff – the function of Hoffman,
Ca – Capillary number based on the contact line velocity,
F_hoff^-1 – inverse of the function of Hoffman,
theta0 - advancing or receding contact angle (depending on the advancing or receding motion).

The function must be solved iteratively due to the appearance of the inverse of the function as one of its own arguments. In the implementation, a C++ function object is introduced defining the inverse Hoffman function, which is solved using RiddersRoot (src/ODE/findRoot). ThetaD is then evaluated using Ca number and the value of F_hoff^-1. The Ca number is evaluated using the existing model for contact line velocity from dynamicContactAngle in OpenFOAM.

The interFoam application is slightly modified. The createFields.H containes initialization of mixture viscosity and surface tension coefficient, which are passed as volScalarFields to the dynamicKistlerContactAngle class. The mixture viscosity is updated during time steps. The application with the new dynamicKistlerContactAngle class and a sample damBreak case are here:




and a couple of snapshots from damBreak case (on a still relatively coarse grid):






This model for the dynamic contact angle (but with another simpler model for contact line velocity) has been successively used in other numerical codes in predictions of drop impact. Examples can be found in:

@Article{Sikalo_etal05,
author = {S. Sikalo and H. D. Wilhelm and I. V. Roisman and S. Jakirlic and C. Tropea},
title = {Dynamic contact angle of spreading droplets: experiments and simulations},
journal = {Physics of Fluids},
year = {2005},
volume = {17},
number = {6},
pages = {62103},
}

@Article{Roisman_etal08,
author = {I. V. Roisman and L. Opfer and C. Tropea and M. Raessi and J. Mostaghimi and S. Chandra},
title = {Drop impact onto a dry surface: Role of the dynamic contact angle},
journal = {Colloids and Surfaces A},
year = {2008},
volume = {322},
pages = {183--191},
}

I am aware that this might not be the most elegant or general way of implementation. Therefore any further suggestions for improvement are welcome.

With kind regards,
Edin Berberovic.

[eberberovic]

Sorry, the files were too large for upload.
Here is the interFoam solver with only the part containing the dynamicKistlerGammaContactAngle class. It should be put into transportModels/interfaceProperties/gammaContactAngle/dynamicKistlerGammaContactAngle

interFoamKistler.tgz
damBreak_thetaA_115.tgz

[gerado]

Hi Edin ,
I am also implementing some models of dynamic contact angle in Openfoam . why you didn t use interfaceproperties.C in your model

thanks
gerard

[eberberovic]

Dear Gerard,

I didn't want to create any objects of type interfaceProperties within the class of type ...FvPatchScalarField just in order to get some data from them. In such case, I believe I would have to use also the twoPhaseProperties objects.

I wanted to pass the values for mu and sigma from the application itself (interFoam).
Therefore I used as a basis the implementation of totalPressure BC.

Regards,
Edin.

[gerado]

Dear Edin ,
I thank you for your fast answer . you know in interfaceproperties.C the contact angle ist corrected (correctgammacontactangle)at each time step . where did you correct it in your inplementation wenn you don t use interfaceproperties.

thanks

gerard

[eberberovic]

Dear Gerard,

OK, now I see what you are asking.

Yes, I do use the same function from interfaceProperties.C to correct the contact angle. What I did is the following: I derived the dynamicKistlerGammaContactAngleFvPatchScalarField from gammaContactAngleFvPatchScalarField, the same way as the existing dynamicGammaContactAngleFvPatchScalarField was derived. The same function correctContactAngle from interfaceProperties.c is used. Only the way how thetaD is calculated has changed (I had to put in a function object and make use of RiddersRoot) as well as the way of passing the arguments to the constructor of dynamicKistlerGammaContactAngleFvPatchScalarField.

I did not make it global, I copied the transportModels folder from src into my interFoam application folder, and made an additional folder called dynamicKistlerGammaContactAngle within.

This is what I wanted to post here, but it was to large, so if you want I can email it to you.

Regards,
Edin.

[gerado]

Dear Edin,
I am really happy to have discussed it with you . it would be very nice if you can sent it to me . I have a last question to you . in your code you say
if uwall>0 theta0 = thetaR
if uwall<0 theta0 = thetaA
is that not the opposite?

thanks
Gerard

e-mail mongerard2000@yahoo.fr

[asaha]

Hello Edin,

Thanks for sharing your implementation on dynamic contact angle in interFoam. I have the following question:

Will it be possible to modify the equation for dyanamic contact angle in your implementation with the following equation

|tan(theta)| = a*Ca^(1/3) -
b*(lambda)^0.04*Ca^0.293; a = 7.48, b = 3.28, lambda = 10E-8, Ca = Capillary number

I would be extremely happy to hear from you in this regard.

Thanks again.

A A Saha.

[eberberovic]

Gerard,

No, I believe that this is correct. This is the most general approach I could make using the exsisting model for contact line velocity in OpenFOAM.
uwall (small u) is evaluated as a dot product of nWall and Uwall. nWall is calculated using nHat by subtracting from it the component normal to the wall. So nWall is parallel to the wall and points into liquid (due to normalized grad(gamma) in nHat). Similarly, Uwall is calculated by taking U from the first cell near the wall and subtracting its normal component. So Uwall is also parallel to the wall, but points in the direction of motion.
Now if the liquid advances, these two vectors point in opposite sides, and if the liquid recedes the point to the same side. Therefore the advancing motion is determined by uwall<0>0.
I also put in some lines to print out all these values at the wall. Uncomment them and you will see that this is correct at the wall.

Saha,

I think this should not be a problem. After you calculate the Ca number, you can simply use your expression for theta. In this case you will not need any function objects, since, as far as I can see, this is an ordinary linear equation which you can evaluate directly (i.e. by theta = arctan(...)).
So, remove the function object I put in, and simply replace the expression for theta.

Regards,
Edin.

[eberberovic]

Sorry for the mistype:
Now if the liquid advances, these two vectors point in opposite sides, and if the liquid recedes they point to the same side. Therefore the advancing motion is determined by uwall<0>0.

[eberberovic]

The advancing motion is determined by uwall<0.
The receding motion is determined by uwall>0.

[asaha]

Hello Edin,

I would appreciate if you can share you code and send the same to asaha@pec.edu.

Thanks,

a a saha

[feijooos]

First, this is a nice addition to the OF code. Second, I am trying to get this code to work with OF 1.5. I am getting some errors regarding the -llduSolver. Did anybody successfully implement this?

Thanks for the help!

[farhagim]

Hello Edin,

thanks for sharing your information. I could not open your files. Would you please share your code again or is it possible to help me to implement the new one?

farhangi.mehran@gmail.com

Originally Posted by eberberovic

The advancing motion is determined by uwall<0.
The receding motion is determined by uwall>0.

[benni]

I try to implement Kistler's dynamic contact angle model in OpenFOAM 2.3.0. With a few changes in Edin's files, I almost could get it work.
Can anybody send me his implementation and particularly the files of a case?
uaecf@student.kit.edu

Thanks a lot!

[yaoyina]

Hello Edin,

thanks for sharing your information. I could not use your files to implement the Kistler model in OF-2.3.0. Would you please share your code again or is it possible to help me to implement the new one?

Thank you very much!

yao-yn14@mails.tsinghua.edu.cn

[kmou]

Hi,

this is all very interesting. Thank you for sharing with us. I am having the same issue as many others, unable to open the .unk folders.
I would greatly appreciate it if you could send it via email, camille.bilger@free.fr
Thank you

[kmou]

Hi Benni,
Do you have the original files of Edin?

[ageorg]

I could open the files normally using archive manager in ubuntu. Now I have the same problem I have to translate the implementation in a newer version of OpenFOAM? I will try the following days and keep you informed. In the meanwhile did anyone managed such a version translation?

Thank you very much in advance

A.

[kmou]

Originally Posted by ageorg

I could open the files normally using archive manager in ubuntu. Now I have the same problem I have to translate the implementation in a newer version of OpenFOAM? I will try the following days and keep you informed. In the meanwhile did anyone managed such a version translation?

Thank you very much in advance

A.

Yes I did put it in OpenFoam 2.3, it works - as in compiles and runs and prints to the screen the dynamic contact angle, but I did it quickly and "blindly" so I still need to understand the interface between the solver and the B.C. to make sure all the variables are passed correctly etc.