View Issue Details

IDProjectCategoryView StatusLast Update
0002864OpenFOAMFeaturepublic2018-10-19 01:18
ReporterSamMallinsonAssigned Tohenry 
PrioritynormalSeverityminorReproducibilityalways
Status resolvedResolutionfixed 
PlatformGNU/LinuxOSDebianOS Versionbuster
Product Version5.0 
Fixed in Versiondev 
Summary0002864: Request for change in direction convension in dynamicAlphaContactAngle BC
DescriptionRegarding the dynamicAlphaContactAngle boundary condition, I think there are a couple of issues:
1. the limiting behaviour of the dynamic contact angle (DCA) does not seem right:

DCA = theta0_ + (thetaA_ - thetaR_)*tanh(uwall/uTheta_)

for uwall >> 0, DCA = theta0_ + (thetaA_ - thetaR_)
for uwall << 0, DCA = theta0_ - (thetaA_ - thetaR_)

2. the value of uwall is such that uwall > 0 is for flow into the liquid, ie DCA = thetaR_; and for uwall < 0, DCA = thetaA_. This means that when the liquid is advancing, it will adopt the receding contact angle, whilst when it is receding, it will adopt the advancing contact angle. This could be fixed, notwithstanding concern #1, by changing the sign of uwall, either in the expression where it is calculated, or by changing the sign of the arguement to tanh.
Steps To ReproduceRun a droplet impacting a surface case: here I have used a 2D case to speed up the runs, but I have also run in 3D and got the same behaviour. Case tar-balls are attached.

1. constantAlphaContactAngle with theta = 120 deg (advancing.tar)
2. constantAlphaContactAngle with theta = 45 deg (receding.tar)
3. dynamicAlphaContactAngle with thetaA = 120 deg, thetaR = 45 deg and thetaE = 90 deg (dynamic.tar)

(I chose these as typical values, see https://www.researchgate.net/profile/Ilia_Roisman/publication/264855677_7_th_Impacting_Droplets_Dynamic_Contact_Angle_Modeling_in_OpenFOAM_R/links/545107b50cf285a067c67f0f/7-th-Impacting-Droplets-Dynamic-Contact-Angle-Modeling-in-OpenFOAM-R.pdf)
Additional InformationAn animation comparing the results from the three cases is also attached. You can see how for the dynamic case, as the droplet spreads, it adopts a low contact angle (receding) and whilst the droplet retracts, it adopts a high contact angle (advancing); note also the the advancing and receding angles are not those specified, but instead as per Description.
TagsNo tags attached.

Activities

SamMallinson

2018-03-01 05:54

reporter  

advancing.tar (30,720 bytes)

SamMallinson

2018-03-01 05:54

reporter  

dynamic.tar (30,720 bytes)

SamMallinson

2018-03-01 05:55

reporter  

receding.tar (30,720 bytes)

SamMallinson

2018-03-01 05:59

reporter  

droplet2d_adv_vs_rec_vs_dyn.gif (3,302,915 bytes)

SamMallinson

2018-03-01 06:01

reporter   ~0009356

see also discussion at cfd-online:

https://www.cfd-online.com/Forums/openfoam-solving/58337-openfoam141dev-new-implementation-dynamickistlergammacontactangle.html#post188397

plus two subsequent posts

henry

2018-03-01 14:25

manager   ~0009364

The definition of advancing and receding depends on the phase these are defined for, if defined for the other phase they would need to be reversed.

SamMallinson

2018-03-04 23:46

reporter   ~0009368

I agree that you need to be consistent with definitions. I used the contact angle definition implicit in the OpenFOAM-5.x/interFoam tutorial capillaryRise, with theta being defined such that it represents the contact angle of the liquid phase with the wall.

When you view the animated .gif attached to this issue, you will see that the dynamicAlphaContactAngle BC mimics the receding angle (45 deg) as the liquid spreads out (advances), and the advancing angle (120 deg) as the liquid retracts (recedes).

With regards the first issue (asymptotic behaviour), I've attached a spreadsheet demonstrating my concern.

dACA.xlsx (15,741 bytes)

henry

2018-03-05 09:59

manager   ~0009369

Currently the definition of the interface direction is with respect to the interface normal:

    // Calculate Uwall resolved normal to the interface parallel to
    // the interface
    scalarField uwall(nWall & Uwall);

where nWall it obtained from the interface normal which in turn is obtained from the interface gradient. This is all internally consistent.

My understanding as that you wish to define the interface normal to be in the opposite direction which I also understand but this change would cause all existing cases using this BC to run incorrectly. How would you propose to support existing cases or would you require all current users of this BC to change all existing cases to conform to your preferred convention?

SamMallinson

2018-03-06 04:34

reporter   ~0009372

Maybe I misunderstand the convention for contact angle definition in dynamicAlphaContactAngle. Are you saying that the convention is different to constantAlphaContactAngle? That is:

- if I wish to specify an advancing contact angle of 120 deg, as defined according to the convention for constantAlphaContactAngle, I should use a value 180 - 120 = 60 deg in dynamicAlphaContactAngle?

- similarly, for a receding angle of 45 deg, as defined according to the convention for constantAlphaContactAngle, I should use a value 180 - 45 = 135 deg?

If this is the case, problem solved.

Perhaps the .H file for this BC could include this detail to avoid confusion.

BTW, the attached animated .gif shows what happens when I assume that the specification of contact angle in dynamicAlphaContactAngle follows the same convention as per constantAlphaContactAngle.

henry

2018-03-06 08:02

manager   ~0009373

Currently the definition of the interface direction is with respect to the interface normal:

    // Calculate Uwall resolved normal to the interface parallel to
    // the interface
    scalarField uwall(nWall & Uwall);

where nWall it obtained from the interface normal which in turn is obtained from the interface gradient. As you see in the above code this normal is dotted with the velocity and when the result is positive the advancing contact angle is used and when negative the receding contact angle is used. The concept of advancing and receding is not present in the constantAlphaContactAngle.

I don't mind changing the convention to your proposal and see your argument for it, however, all existing cases would no longer work correctly and I am asking you how this should be handled.

SamMallinson

2018-03-08 05:07

reporter   ~0009390

> all existing cases would no longer work correctly

I suspect that most users would expect, as I did initially, to define contact angles for the dynamicAlphaContactAngle BC in the same way as they are defined for the constantAlphaContactAngle BC.

That is, the existing cases were likely not properly specified: users would most likely have defined an advancing contact angle of 120 deg by setting thetaA as 120, not the supplementary angle, 60 deg, etc.

I tried to highlight this using the examples and the animation.

So, I would suggest you make the definitions consistent, thus 'breaking' (or more likely, fixing) the existing cases.

I did this (see attached .tar file - I called the modified BC dynamicAlpha1ContactAngle) by changing the calculation of theta0_ from:

return theta0_ + (thetaA_ - thetaR_)*tanh(uwall/uTheta_);

to:

return theta0_ - (thetaA_ - thetaR_)*tanh(uwall/uTheta_);

dA1CA.tar (20,480 bytes)

SamMallinson

2018-03-08 06:00

reporter   ~0009391

See also attached animated .gif which compares the current method (dynamicAlphaContactAngle) with the proposed method (dynamicAlpha1ContactAngle).

dACA_original_vs_proposed.gif (3,731,459 bytes)

henry

2018-03-08 07:14

manager   ~0009392

All the cases I have are properly defined and work correctly with the current implementation. With your proposed change in direction convention they would still run but the results would be incorrect. This would also be the case for anyone else who had setup cases correctly with the current convention.

guin

2018-03-08 10:59

reporter   ~0009394

Just to add the point of view of a third one (mine). As a user, I agree that the current convention is confusing and might lead to wrong setup of new cases. In general, I would expect advance contact angles to be bigger than receding ones.

https://en.wikipedia.org/wiki/Contact_angle#/media/File:Contact_angle.svg
https://en.wikipedia.org/wiki/Contact_angle#/media/File:Dynamic_contact_angle_measurement.svg

If the main for not accomplishing changes is just to avoid breaking compatibility with prior setups, then we should evaluate the priorities (i.e. existing cases vs newcoming ones)... I would say that people that previously set up simulations generally have more experience than people doing it for the first time. And if this change goes towards a more intuitive formulation, they may even be glad of that.

henry

2018-03-08 13:01

manager   ~0009395

I have used this BC I wrote in the late 90's successfully for many years. It has also been used by many people over that period and I do not think in reasonable to propose a change which would break all existing cases or suggest that they do not know how to use this BC and that their cases are probably incorrect anyway.

I think that any change to this BC should provide for backward-compatibility.

SamMallinson

2018-03-12 00:00

reporter   ~0009406

@guin - thanks for your comment

@henry - I take your point. How about an extra parameter (I'll call it "physical" but maybe you can think of a better name), which switches between 0, if one wishes to use your original defintion, and 1, if one wishes to use the physical definition, of contact angle.

The conditional test to differentiate between these two values is:

    scalar supplement_ = 1.0; // use supplementary angle definition
    if (physical_ == 1.0) // otherwise use physical angle definition
    {
        supplement_ = 0.0;
    }

and the return call is:

    return theta0_ + (supplement_ - physical_)* (thetaA_ - thetaR_)*tanh(uwall/uTheta_);

I've attached the BC code (I've left it called dynamicAlpha1ContactAngle, please change it to dynamicAlphaContactAngle if your choose to adopt it).

dA1CA_proposed.tar (20,480 bytes)

SamMallinson

2018-03-12 00:03

reporter   ~0009407

Attaching test case for which physical = 0 (ie using original definition of contact angle in dynamicAlphaContactAngle).

Please feel free to use this as a tutorial case. You may also like to change the definition of values in the file system/parameters, or in 0.org/alpha.water

Also, there is a vestigial codeStream in system/controlDict for calculating the kinetic energy of the flow - this can be removed if you would prefer.

supplementary.tar (30,720 bytes)

henry

2018-03-12 00:29

manager   ~0009408

Is the switch you propose a scalar as indicated in

if (physical_ == 1.0)

if so why? Is a value other than 1.0 or 0.0 valid? If not shouldn't it be a boolean?

SamMallinson

2018-03-12 00:36

reporter   ~0009409

I noticed that I had not given the correct definition of physical in one of the constructors: I've changed it to

physical_(dict.lookupOrDefault<scalar>("physical", 0))

Attaching correct (ie including the above) proposed BC.

Now the case will run without physical being defined in the dict for the BC.

dA1CA_proposed_correct.tar (20,480 bytes)

SamMallinson

2018-03-12 00:48

reporter   ~0009410

@henry - thanks for the prompt response.

Of course, it should be boolean, my mistake. I've changed the expression in the constructor to:

    physical_(dict.lookupOrDefault<bool>("physical", 0))

and the conditional statement to:

    bool supplement_ = 1; // use supplementary angle definition
    if (physical_ == 1) // otherwise use physical angle definition
    {
        supplement_ = 0;
    }

SamMallinson

2018-03-12 01:27

reporter   ~0009411

Attaching animation comparing results for physical = 0 (ie original CA definition) and physical = 1 (ie constantAlphaContactAngle CA definition).

dA1CA_supp_vs_phys.gif (4,036,708 bytes)

SamMallinson

2018-03-12 01:28

reporter   ~0009412

Attaching BC code which had physical as a boolean parameter, dA1CA_boolPhys.tar.

dA1CA_boolPhys.tar (20,480 bytes)

henry

2018-03-12 09:04

manager   ~0009413

I don't agree with your contention that your proposed convention is "physical" implying the current convention is not; it is just a convention.

I think the choice of convention should be handled in the input, there is no need for a complex switch and blending statement in the implementation.

I will write a simpler version which gives the user a choice of input convention without adding clutter to the actual implementation of the BC.

SamMallinson

2018-03-12 22:03

reporter   ~0009415

Hi @henry,

> I don't agree with your contention that your proposed convention is "physical" implying the current convention is not; it is just a convention.

Perhaps "consistent" as in consistent with other conventions in OF? Anyway, it is just a name, and you'll be re-writing the code.

My main intention was to make the change transparent to earlier users of this BC, hence having a switch which did not need to be present for the modified code to work for the earlier users. This seemed to be your main concern in our earlier discussions.

> I think the choice of convention should be handled in the input

Yes, that's what was done, with the use of an extra switching parameter, in the BC dict.

> I will write a simpler version which gives the user a choice of input convention without adding clutter to the actual implementation of the BC

I tried my hardest to make the change minimal, and make it as simple as possible: one extra parameter, one extra if statement and a minor change to the return statement. I'm looking forward to seeing how your code looks.

Cheers,

Sam

gdmcbain

2018-03-13 03:41

reporter   ~0009417

A suggestion: In the forthcoming "simpler version", perhaps the two parameters thetaA_ and thetaR_ which appear in the implementation only via their difference could be replaced by a single parameter, say deltaTheta; i.e. replace

    return theta0_ + (thetaA_ - thetaR_)*tanh(uwall/uTheta_);

with

    return theta0_ + deltaTheta*tanh(uwall/uTheta_);

Then the argument about which phase was advancing or retreating would be moot, being replaced by an arbitrary convention on the sign of the new parameter. I admit that initially I too had thought that the advancing phase was that same one for which the contact angle had been defined.

Legacy examples could be supported by an alternate constructor which formed deltaTheta = thetaA_ - thetaR_.

henry

2018-07-17 14:06

manager   ~0009859

Resolved by commit a9cdfa0f4e953a11839502c284501c668baea231

chris

2018-07-17 15:15

manager   ~0009861

OpenFOAM must be funded:
https://openfoam.org/news/funding-2018/

This kind of prioritisation over features reported through the issue tracking system requires a silver support package:
https://openfoam.org/maintenance/

The cost is a small fraction of a single, annual, commercial CFD licence. It is a small fraction of the total cost of any serious, in-house CFD capability.

Issue History

Date Modified Username Field Change
2018-03-01 05:54 SamMallinson New Issue
2018-03-01 05:54 SamMallinson File Added: advancing.tar
2018-03-01 05:54 SamMallinson File Added: dynamic.tar
2018-03-01 05:55 SamMallinson File Added: receding.tar
2018-03-01 05:59 SamMallinson File Added: droplet2d_adv_vs_rec_vs_dyn.gif
2018-03-01 06:01 SamMallinson Note Added: 0009356
2018-03-01 14:25 henry Note Added: 0009364
2018-03-04 23:46 SamMallinson File Added: dACA.xlsx
2018-03-04 23:46 SamMallinson Note Added: 0009368
2018-03-05 09:59 henry Note Added: 0009369
2018-03-06 04:34 SamMallinson Note Added: 0009372
2018-03-06 08:02 henry Note Added: 0009373
2018-03-08 05:07 SamMallinson File Added: dA1CA.tar
2018-03-08 05:07 SamMallinson Note Added: 0009390
2018-03-08 06:00 SamMallinson File Added: dACA_original_vs_proposed.gif
2018-03-08 06:00 SamMallinson Note Added: 0009391
2018-03-08 07:14 henry Note Added: 0009392
2018-03-08 10:59 guin Note Added: 0009394
2018-03-08 13:01 henry Note Added: 0009395
2018-03-08 22:10 henry Severity major => minor
2018-03-08 22:10 henry Category Bug => Feature
2018-03-08 22:10 henry Summary Improper behaviour of dynamicAlphaContactAngle => Request for change in direction convension in dynamicAlphaContactAngle BC
2018-03-12 00:00 SamMallinson File Added: dA1CA_proposed.tar
2018-03-12 00:00 SamMallinson Note Added: 0009406
2018-03-12 00:03 SamMallinson File Added: supplementary.tar
2018-03-12 00:03 SamMallinson Note Added: 0009407
2018-03-12 00:29 henry Note Added: 0009408
2018-03-12 00:36 SamMallinson File Added: dA1CA_proposed_correct.tar
2018-03-12 00:36 SamMallinson Note Added: 0009409
2018-03-12 00:48 SamMallinson Note Added: 0009410
2018-03-12 01:27 SamMallinson File Added: dA1CA_supp_vs_phys.gif
2018-03-12 01:27 SamMallinson Note Added: 0009411
2018-03-12 01:28 SamMallinson File Added: dA1CA_boolPhys.tar
2018-03-12 01:28 SamMallinson Note Added: 0009412
2018-03-12 09:04 henry Note Added: 0009413
2018-03-12 22:03 SamMallinson Note Added: 0009415
2018-03-13 03:41 gdmcbain Note Added: 0009417
2018-07-17 14:06 henry Assigned To => henry
2018-07-17 14:06 henry Status new => resolved
2018-07-17 14:06 henry Resolution open => fixed
2018-07-17 14:06 henry Fixed in Version => dev
2018-07-17 14:06 henry Note Added: 0009859
2018-07-17 15:15 chris Note Added: 0009861