open TELEMAC-MASCARET The mathematically superior suite of solvers

In respect to the topic Attached please find example with modified DRAGFO subroutine that allows a number of polynomial drag force zones to be defined and the dragforce stress applied to these zones to be calculated for both linear and quasi-bubble discretization

Regards
Tony

Hi Tony & JMH,

Thanks for your replies, however I am still finding discrepancies in my results. To demonstrate this I have attached the run files for 2 simulations of flow through a simple channel. Both use exactly the same geometry and boundary conditions files, the only difference between the two is the dragfo code. In the first no calculation is made for AIRE, in the second, the calculation for AIRE is included. I have corrected the write UNORM statement so that U.txt now opens and closes outside of the do loop based on JMH's instruction (however this made no difference to how the results were written to U.txt). On a quick inspection of the results in U.txt, I see that the max velocity using dragfo_no_aire.f90 is 1.8695m/s, whilst the maximum velocity using dragfo_with_aire is 1.8669m/s. Whilst the difference is small, why is there a difference at all??? I have other eamples where the difference is greater but have used this example for its simplisity.
Any insight would be greatly appreciated.
Cheers
Danny

Hello,

If realise now that if you look at the original dragfo.f subroutine, you have in the end:

IF(AIRE.GT.1.D-6) THEN
SOM = 1.D0 / AIRE
ELSE
IF(LNG.EQ.1) WRITE(LU,*) 'DRAGFO : AIRE DE LA ZONE NULLE'
IF(LNG.EQ.2) WRITE(LU,*) 'DRAGFO: AREA OF ZONE EQUAL TO ZERO'
CALL PLANTE(1)
STOP
ENDIF
!
CALL OS('X=CX ',X=FUDRAG,C=SOM)
CALL OS('X=CX ',X=FVDRAG,C=SOM)

which is necessary to get the correct value of FUDRAG and FVDRAG, this is why we compute the area AIRE and I do not understand where this is in your implementation.

Then if you look at your Fortran files with AIRE or without AIRE, in one there is the section for quasi-bubble velocity, in the other it is not there. As you have a quasi-bubble velocity in your case there is no wonder that this makes a difference. You should have the whole stuff, section for quasi-bubble followed by the division by AIRE.

With best regards,

Jean-Michel Hervouet

Hi Jean-Michel,

Thanks for your reply, this makes sense. Just to confirm, does the attached dragfo look ok to you?

Cheers

Danny

Hello,

Well I do not really understand what is your theory, why multiplying by AIRE and then dividing by AIRE, why two lines with AIRE=0.D0, what is the use of AIRE_B, why the division by the depth ?

The basic subroutine gives you the example of one pier with diameter DIAM, in the formula of the drag force you have an area that should be DIAM*depth, as in the non conservative equations there is a division by the depth this leaves only DIAM, so the formula of the drag force gives -0.5*DIAM*CD*UNORM, this force is spread over an area AIRE as a local stress, so the division by AIRE.

With best regards,

JMH

Hi JMH,
I apologies, I realise that wasn't the clearest I could have been. To better define my problem let me start again using the example in your book for an array of turbines, I want to express the drag from an array of turbines in terms of the array density k, where k= total swept area/array plot area.

As a starting point, the array force is given as:

F=0.5 * rho * U^2 * Cd * A * n

Where A is the swept area of a single turbine and n is the number of turbines in the array. To express as a shear, divide through by AIRE:

F=(0.5* rho * U^2 * Cd * A * n) / (AIRE)

Now substitute in array density, k, where k=(A*n)/AIRE:

F=0.5* rho * U^2 * Cd * k

Since the turbines do not take up the whole depth, and then divide through by depth, divide through by rho and UNORM to give the final expression for FUDRAG and FVDRAG:

F=(0.5 * U * Cd * k) / depth -EQ A

By taking this approach I need to order the operations in my DRAGFO file as follows (as shown also in attachment):

1. Loop to calculate AIRE
2. Correct AIRE for qausi bubble discretization
3. to calculate FUDRAG using EQ A

As such there is no need to include SOM because the array density k incorporates 1/AIRE. Does this approach look ok to you?



Cheers
Danny

Hello,

OK, that makes sense now. Because the turbines do not occupy the whole depth, you fall back on a formula that looks more like a friction law, where you have the division by the depth because we deal with the non conservative momentum equations. This means also that a specific local friction coefficient would do an equivalent job.

With best regards,

Jean-Michel Hervouet

Hi Jean-Michel,

Yes that was my previous approach. So to confirm does the dragfo file (Dragfo_telemac_forum.txt) look ok to you?

Cheers

Danny

Yep !

JMH

Hi Pham!

I am trying to use drag force subroutine to simulate a wall. I get the following error after trying to run the cas.
DRAGFO: AREA OF ZONE EQUAL TO STOP 1

what is the meaning of the error? I tried to fix the polygon points and everything looks okay. I am wondering where to fix it now! I attached the screenshot. Do you have any thoughts?