Learning Forum

[nahian.masud]

Is there any way to find the lift force in CG motion udf? I need to to use in the vertical velocity. Here is the code I have written so far. The 2D airfoil is flapping, so I described a flapping motion. The airfoil has boundary ID of 7. now I use Compute_Force_and_Moment macro to find lift force and use it in the vertical velocity. 
The problem is, it says floating point error after 2 or 3 time steps ( residuals goes super high and simulation stops), but if i dont use the compute force and moment macro, the flapping only works well. I need to find the lift force anyhow to get the transational velocity. can anyone help me please?

#include "udf.h"

DEFINE_CG_MOTION(paperfreq,dt,vel,omega,time,dtime)

{

Thread *t; 

face_t f;

/* get moment or force */

Domain *d = Get_Domain(1); /* for single phase flow*/

 

Thread *t_object = Lookup_Thread(d, 7); /* you get Boundary_ID from the boundary condition panel in Fluent*/

real force[ND_ND], moment[ND_ND], cg[ND_ND]; /*initialise*/

 

Compute_Force_And_Moment (d,t_object,cg,force,moment,TRUE);

 

 /* force components of surface "Boundary_ID" real force_y = force[1];*/

 

real w,a, pi; 

pi = 3.14159265; 

a = (pi*15)/180; /*pitch amplitude*/ 

w = 0.62832; /*angular velocity*/  

/*defining the flapping motion*/ 

omega[2]=-a*2*pi*2.97*cos(2*pi*2.97*time); /*pitching motion*/

vel[1]= force[1]*dtime/0.0116;

 

}

The code above shows the error with very high residuals. But the bottom code works just fine,

#include "udf.h"

DEFINE_CG_MOTION(paperfreq,dt,vel,omega,time,dtime)

{

Thread *t; 

face_t f;

real w,a, pi; 

pi = 3.14159265; 

a = (pi*15)/180; /*pitch amplitude*/ 

w = 0.62832; /*angular velocity*/  

/*defining the flapping motion*/ 

omega[2]=-a*2*pi*2.97*cos(2*pi*2.97*time); /*pitching motion*/

 

}

UPDATE:
I found the possible reason of the error (not sure yet)

Thread *t_object = Lookup_Thread(d, 7); /* Boundary_ID from the boundary condition panel in Fluent*/

In this line the ID of my airfoil is 7. It is a wall. If I use this, the error occurs. but if I use other ID, the problem does not occur. But I need to use the airfoil in order to find the lift force on it.

[Rob]

If you report the force on the airfoil wall (print the return for force[1] ) for some timesteps but use the second block of commands in the UDF what sort of numbers are you seeing relative to what you expect? Also, I can't see a mass in the system - what stops the wing just accelerating upwards? 

[nahian.masud]

Hi, thanks for replying. this is the latest UDF I am using,

#include "udf.h"

#define MASS 0.0143

 

/* Declare a global variable to store the lift force */

static real lift = 0.0;

 

/* Function to compute the CG motion, including flapping motion */

DEFINE_CG_MOTION(paperfreq, dt, vel, omega, time, dtime)

{

    Domain *d = Get_Domain(1); // Assuming single phase flow

    Thread *t_object = Lookup_Thread(d, 7); // Airfoil's thread ID

    real force[2], moment[2], cg[2]; // Initialize arrays

    real x=0;

    NV_S(vel, =, 0.0);

    /* Compute force and moment at the boundary */

    Compute_Force_And_Moment(d, t_object, cg, force, moment, FALSE);

 

    /* Update global variable for lift force */

    lift = force[1];

 

    real w, a, pi;

    pi = 3.14159265;

    a = (pi * 15) / 180; // Pitch amplitude in radians

    w = 0.62832; // Angular velocity in rad/s

 

    /* Define the flapping motion */

    omega[2] = -a * 2 * pi * 2.97 * cos(2 * pi * 2.97 * time); // Pitching motion

    vel[1]=lift*dtime/MASS;

    x=vel[1];

    

Message("lift:%f x:%f\n",lift,x);

}

The problem is, if i use the vel[1], the lift force goes very high, so does the velocity. This makes the airfoil go outside domain and it creates floating point exception error. If I dont use the vel[1], and I print the lift force using Message, it shows correct lift value. 

What I need is, change the location of the airfoil based on the lift force. If there is any other way to give it a y-displacement it is also ok. I am giving it a flapping motion, which generates various lift force based on time and I need the latest lift force to update the airfoil location before going to next time step.

[Rob]

OK, so the reports are correct, the macros are doing what you expect but you're moving the airfoil too far with the force? Would reducing the speed or time step help? 

[nahian.masud]

 

yea I reduced the time step. Usually the the velocity range should be below 0.32 m/s and I made the timestep in a way that, the velocity stays in range. But still the lift force is going super high within a few timesteps. I cannot find the reason for it.
and yes, I am moving the airfoil too far away with the force, as it looks like.


alternatively, I was thinking that, I can write the force data in .OUT or .txt file using report definition. Is there any way, to read the force data from that .OUT file and use it in the CG_motion?

thanks

 

[Rob]

How far is the airfoil moving relative to the cell size in one timestep? 

[nahian.masud]

This is what I am getting.
Is the following way possible?I can write the force data in .OUT or .txt file using report definition. Is there any way, to read the force data from that .OUT file and use it in the CG_motion for every time step?

[Rob]

You'd use the value directly into the CG_MOTION macro, no need to use a file. 

[nahian.masud]

There is no single constant value. The force is changing every timestep due to changing angle of attack (as it has oscillating pitching motion). So how can I use the updated force value, that is generated from the fluent solver?

Also can you please tell me, how to modify the vel[1] so that my lift force does not go super high?

thanks

[Rob]

The CG_MOTION examples should cover the data use - basically the force will be calculated and then the motion assigned. 

Regarding vel[1] values. At present it's uncertain if they're high because the maths you're using isn't quite right or if it's because the solver is diverging and that then returns a nonsense value to the equation which in turn causes the airfoil to move too far.