[dtmith]

Hi all,

I've been running LES simulations of plate rotating about an axis. It starts from rest, accelerates to a certain velocity and then continues to rotate at that velocity for the rest of the simulation.

I define an "omega" profile for the acceleration portion. It is made up of a quartic polynomial followed by a linear section followed by another quartic. I noticed in the force data a non-physical spike when it transitioned from the quartic profile to the linear profile and back again.

Here's my UDF:

#include "udf.h"

DEFINE_ZONE_MOTION(motion3,omega,axis,origin,velocity,time,dtime)
{
	real pi, w, tscl;
	
	pi = 3.141592654;
	w = -1.402;
	tscl = 5.58505360638183;
	
	if (time<0.1/tscl)
	{
		*omega = (-5000/9*pow((time*tscl),4.) + 1000/9*pow((time*tscl),3.))*w;
	}
	else
	{
		if (time<0.9/tscl)
		{
			*omega = (10/9*(time*tscl) - 1/18)*w;
		}
		else
			if (time<1/tscl)
			{
				*omega = (-(-5000/9*pow((1-(time*tscl)),4.) + 1000/9*pow((1-(time*tscl)),3.))+1)*w;
			}
			else
			{
				*omega = w;
			}
			
	}
	Message ("time = %f, omega = %f\n", time, *omega);
	N3V_D(axis,=,0.0,1.0,0.0);

	}

So basically w and tscl are scaling factors for omega and time respectively. For time is 0 to 0.1/tscl, my velocity profile is this quartic...

*omega = (-5000/9*pow((time*tscl),4.) + 1000/9*pow((time*tscl),3.))*w;

Then from 0.1/tscl to 0.9/tscl the profile is linear...

*omega = (10/9*(time*tscl) - 1/18)*w;

Then from 0.9/tscl to 1/tscl it goes to another quartic...

*omega = (-(-5000/9*pow((1-(time*tscl)),4.) + 1000/9*pow((1-(time*tscl)),3.))+1)*w;

Finally that quartic smoothly transitions to a constant rotation velocity above 1/tscl

*omega = w;


After a LOT of screwing around I discovered the problem was with the linear portion, typing my constants as fractions (10/9 and 1/18) meant the fractions were ignored, it was basically evaluating...

(time*tscl)*w instead of (10/9*(time*tscl) - 1/18)*w

If I type out the fractions as 1.1111111 and 0.05555555 instead of 10/9 and 1/18 it seems to work fine.

But for the other segments, I am also using fractions to define omega (5000/9) and it seems to be working fine?

Does anyone know why for that particular line it can't accept a fraction? I've looked at it so many times and can't figure out why it would work for the quartic but not the linear portion. Maybe it'll save me from screwing up in future, lol.

[`e`]

Well for a start, add a trailing period to your integer numerators and denominators. For example: "10./9.". This syntax tells the C compiler to treat these numbers as real data types rather than integers because otherwise the compiler will evaluate this fraction and keep the type as an integer (perhaps 10/9 == 1). The precise treatment of integer divisions depend on the compiler you're using (C89 behaves differently to C99 etc).