## Calculation of forces and moments

Questions and remarks about code_saturne usage
Forum rules
Luciano Garelli
Posts: 245
Joined: Fri Dec 04, 2015 1:42 pm

### Calculation of forces and moments

Hi everybody,

I'm solving a problem of a rotating object (cylinder) in a fluid, initially at rest, and I want to compute the forces and moments acting on it. The object is placed at the origen and a tangencial (sliding) velocity is prescribed at the surface through cs_user_boundary_conditions.f90.

The forces and moments (with respect to the origen) are computed using cs_user_extra_operations-global_efforts.f90, for the giver boundary

Code: Select all

``````! Face center
xnod(1) = cdgfbo(1, ifac)
xnod(2) = cdgfbo(2, ifac)
xnod(3) = cdgfbo(3, ifac)

! Global Efforts
do ii = 1, ndim
xfor(ii) = xfor(ii) + bfprp_for(ii, ifac)
enddo

!Global Moments
tfor(1)= tfor(1)+ (bfprp_for(3, ifac)* xnod(2) - bfprp_for(2, ifac)*xnod(3))
tfor(2)= tfor(2)+ (bfprp_for(1, ifac)* xnod(3) - bfprp_for(3, ifac)*xnod(1))
tfor(3)= tfor(3)+ (bfprp_for(2, ifac)* xnod(1) - bfprp_for(1, ifac)*xnod(2))
``````
As a first test in order to verify the moments computation, the fluid is at rest, so the direction of the resulting moment has to be opposed to the direction to the imposed rotation, isn't it?.

The problem is that the resulting moment direction is not aligned with the rotation direction. Any idea why?

I'm not using any turbulence model and the convergence is ok.

I attach the file where is setted the tangencial velocity and the forces and moments computed.

Attachments
cs_user_extra_operations-global_efforts.f90
cs_user_boundary_conditions.f90
Yvan Fournier
Posts: 3188
Joined: Mon Feb 20, 2012 3:25 pm

### Re: Calculation of forces and moments

Hello,

Did you visualize the stresses/efforts vector fields on the boundary using ParaView, EnSight, or equivalent (just to check that values are realistic/what you would expect) ?

As in incompressible computations, pressure is defined relative to an arbitary constant, you might have tricky to interpret results if the surface you integrate on is not closed. If is is closed, the constant may balance out, so things may be more surprising, but in this case, visualization is the first step to understanding what is not working as expected (or pinpointing a bug if such is suspected).

Also, it is not clear relative to which axis you are computing a moment: it seems you are computing the moment contribution relative to each face's center of gravity, rather than to a reference axis (or I'm missing the formula).

Regards,

Yvan
Luciano Garelli
Posts: 245
Joined: Fri Dec 04, 2015 1:42 pm

### Re: Calculation of forces and moments

Hello Yvan and thanks for answering.

The first step was check the if the velocity field it is what I expected and then the efforts. The velocity was ok, but the efforts looks strange, so I decide to simplify the problem and solve a rotating sphere with the fluid at rest.

For this case the inner sphere has radius r_inner=0.1 [m] and an imposed angular velocity W = (10 i; 0 j; 0 k) [rad/s], wich is set thorugh the cs_user_boundary_conditions.f90 where the velocity at the surface of the sphere is calculated doing V= W x R, being R the vector to the boundary face center. I would expect a maxium velocity of 1 [m/s] and decreasing to the axis of rotation (x axis).

The resulting efforts will be purely viscous and it has to be proportional to the velocity. To compute the moments I did M = R x F, where F are the efforts at the face center. This is implemented in cs_user_extra_operations-global_efforts.f90.

I attach a pair of figures with the velocities, velocity vectors and efforts. The velocity magnitude and directions seems to be ok, but the resulting efforts are quite strange. May be I am setting wrong the boundary condition at the sphere surface?

Regards
Attachments
cs_user_boundary_conditions.f90
Yvan Fournier
Posts: 3188
Joined: Mon Feb 20, 2012 3:25 pm

### Re: Calculation of forces and moments

Hello,

If your case is small enough, could you post the mesh and the setup also ?

I'll try to take a look sometimes next week, to check for bugs (won't have enough time before that).

Regards,

Yvan
Luciano Garelli
Posts: 245
Joined: Fri Dec 04, 2015 1:42 pm

### Re: Calculation of forces and moments

Hello,

The case is small, so I uploaded the mesh, the xml and the subroutines.

In order to find where is the problem (I think in my subroutines) I have solved the same problem but now the inner sphere is held fixed and the sourrounding fluid is initialized with a rotating velocity throught the GUI

Code: Select all

``````mag_omega=10;
omegar=mag_omega;
omegar=0;
omegar=0;
velocity = omegar*z-omegar*y;
velocity = -(omegar*z-omegar*x);
velocity = omegar*y-omegar*x;
``````
In this case the efforts are purely tangencial But in the previous case, with the rotating sphere, the velocity field over the sphere is ok, but not the efforts and appears nomals efforts also. Regards.
Attachments
Sphere_case.tar.gz
Yvan Fournier
Posts: 3188
Joined: Mon Feb 20, 2012 3:25 pm

### Re: Calculation of forces and moments

Hello,

I ran your case for a few time steps, and indeed, the results seem strange (I would not be surprised by difficult-to-interpret constraints due to the pressure component, but I would expect the tangential component to be aligned (or almost) with the velocity.

I'll need to check this next week with colleagues who know this part of the code better than I do...

I'll keep you updated.

Regards,

Yvan
Yvan Fournier
Posts: 3188
Joined: Mon Feb 20, 2012 3:25 pm

### Re: Calculation of forces and moments

Hello,

I think I have found the bug, which is due to an incorrect parameter in 1 the z direction.

Could you try using the attached file (either in code sources in src/base, or in your user subroutines) ?

I'll push the fix on Monday, so it will be included in all future releases...

Regards,

Yvan
Attachments
condli.f90
Luciano Garelli
Posts: 245
Joined: Fri Dec 04, 2015 1:42 pm

### Re: Calculation of forces and moments

Thanks Yvan,

Now the results looks good and the values are as expected.

One additional comment/suggestion for future versions, when choosing a slinding wall in the Boundary Conditions options through the GUI, only can be applied a traslational velocity.

A rotational velocity will be very helpfull and with the definition of the angular velocity [rad/s] and the origin [x0, y0, z0] and direction [xi, yi, zi]of the rotation, the velocity at the wall can be imposed.