Yes maybe using the number of particles on the wall with impingement location can be fed into a UDM for post processing, I'll keep you posted on this.
But before I can go that far, I see that my Force values and pressure values seem to be off. I printed the variable values for 'Fa, Fn, radius, P_local' from my UDF into a text file.
This is how the values look in the same chronological order as stated above (Fa,     Fn,      radius,   P_local) 1.000000 0.023929 0.021373 0.000000 1.000000 0.027982 0.023274 0.000000 1.000000 0.024288 0.021810 0.000000 1.000000 0.027118 0.023750 0.000000 1.000000 0.026723 0.024311 0.000000 1.000000 0.024109 0.021391 0.000000 1.000000 0.026214 0.024212 0.000000 1.000000 0.027531 0.023820 0.000000 1.000000 0.028039 0.023649 0.000000 1.000000 0.025427 0.023606 0.000000 1.000000 0.024226 0.021913 0.000000 1.000000 0.025716 0.022863 0.000000 1.000000 0.025270 0.023653 0.000000.....................& so on
What could be the reason for such random values being calculated? most importantly why are the values of Pressure equal to zero(found from F_P(f, t) = 0)?
I know this is a lot to ask for on a forum like this but your thoughts on this would help me a lot.
Â
UDF:
Â
#include "udf.h"
#include "dpm.h"
#include "math.h"
#include "mem.h"
#include
#define Rho_air 1.225 Â /* Air density in kg/m^3 */
#define Cd 0.47 Â Â Â Â /* Drag coefficient for a sphere */
#define mu 0.5 Â Â Â Â /* Coefficient of friction */
#define Vflow 10
#define PI 3.14
/* Define your boundary condition function */
DEFINE_DPM_BC(particle_bc, p, t, f, normal, dim)
{
  real Vparticle = sqrt(P_VEL(p)[0] * P_VEL(p)[0] + P_VEL(p)[1] * P_VEL(p)[1] + P_VEL(p)[2] * P_VEL(p)[2]);  /* Velocity magnitude in m/s */
  real V = fabs(Vparticle - Vflow);
  FILE *fp1 ;
  FILE *fp2;
  real alpha;  /* angle of particle path with face normal */
  real vn=0.;
  real nor_coeff = 1.;
  real tan_coeff = 0.3;
  //real normal[3];
  int i, idim;
  real NV_VEC(x);
  real radius = P_DIAM(p) / 2.0;
  real A = PI * radius * radius;  /* M_PI is a constant for π */
  /* Get the face pressure */
  real P_local = F_P(f, t);
  // real P_local = C_P(c,t);
  /* Calculate aerodynamic force Fa */
  real Fa = 0.5 * Rho_air * A * Cd * V * V;
  /* Calculate normal force Fn */
  real Fn = P_local * A;
  /* Calculate frictional force Ff */
  real Ff = mu * Fn;
  if (Fa > Ff)
  {
    /* Particle needs to bounce with a coefficient of restitution of 1 */
Â
       /* Compute normal velocity. */
   for(i=0; i
   Â
    vn += P_VEL(p)[i]*normal[i];
   /* Subtract off normal velocity. */
   for(i=0; i
    P_VEL(p)[i] -= vn*normal[i];
   /* Apply tangential coefficient of restitution. */
   for(i=0; i
    P_VEL(p)[i] *= tan_coeff;
   /* Add reflected normal velocity. */
   for(i=0; i
    P_VEL(p)[i] -= nor_coeff*vn*normal[i]; Â
   /* Store new velocity in P_VEL0 of particle */
   for(i=0; i
    P_VEL0(p)[i] = P_VEL(p)[i];
   fp1=fopen ("Impact.txt" , "a") ;
   fprintf( fp1 , "%f %f %f %f %f\n", Fa, Fn, Ff, radius, P_local);
   fclose(fp1);
   return PATH_ACTIVE;
  }
  else
  {
    /* Particle sticks to the boundary */
    for (int i = 0; i < dim; i++)
    {
      P_VEL(p)[i] = 0; /* Set velocity to zero */
    }
    P_USER_REAL(p, 0) = 1.0; /* Optionally set a user-defined variable to indicate sticking */
    /* Save the impingement location */
    P_USER_REAL(p, 1) = P_POS(p)[0];  /* X-coordinate */
    P_USER_REAL(p, 2) = P_POS(p)[1];  /* Y-coordinate */
    P_USER_REAL(p, 3) = P_POS(p)[2];  /* Z-coordinate */
    P_USER_REAL(p, 4) = 1.0;  /* Mark as stuck */
    fp2=fopen ("Impact2.txt" , "a") ;
    // fprintf(fp2 , "%f %f %f %f\n", P_USER_REAL(p, 0), P_USER_REAL(p, 1), P_POS(p)[1], P_POS(p)[2]);
    fprintf( fp2 , "%f %f %f %f %f %f \n", Fa, Fn, Ff, radius, P_local, F_P(f, t));
    fclose(fp2);
    return PATH_ABORT; /* Stop tracking the particle */
  }
}
Â
Â
Â
Thanks, Kiran
Â
Please Login to Report Topic
Please Login to Share Feed
Edit Discussion
You are navigating away from the AIS Discovery experience