[guoxj]

To facilitate the discussion，I open a new theme. Moving the problem description to below:
I am simulating the starting behavior of a domain including Impeller, using a Fortran subroutine.The subroutine is used to update the angular velocity of the domain.And I calculate the angular velocity taking into account the torque generated by the blade parts and the inertia of the impeller, and this value is used by the cfx solver to set the angular velocity of the blade. The subroutine is used to calculate angular velocity of the rotational domain in the timestep n by using the value of previous time step.The problem is, when the calculation is proceeding,although the cycle is set to update angular velocity at the start of a new timestep, the subroutine is called many times per inner iteration,not once every timestep. That is angular velocity has been updated dozens of times within each iteration. I want to what is the reason of it? And how to solve this problem?

The error message shown below:

ERROR #001100279 has occurred in subroutine ErrAction. |
| Message: |
| Angular velocity can depend only on time for transient simulations.

[michael_owen]

The error doesn't have anything to do with you problem description. The angular velocity specified for a rotating domain cannot be a function of anything except time. This is a well known limitation of CEL. You will have to be more crafty about how you do this. Also, FORTRAN is not required.

[michael_owen]

By the way, if you have R12.1 you can do this using the 6DOF rigid body solver. Otherwise you'll have to us technique like the one I developed here:

http://www.youtube.com/watch?v=md8vXAUn7g0

http://www.youtube.com/watch?v=9uAobJm0r1I

These were done in v11, before the 6DOF solver, without FORTRAN.

[guoxj]

Thanks for your reply.But i have not the R12.1,also i can not open the links you give me.could you give me some description about the techniques you have developed in other ways.Or could you give me some other suggestions?

[guoxj]

Here i apology for posting new topic behind Filippo's

[michael_owen]

The technique I developed, and I'm sure I'm not the first, is to use an Additional Variable to keep track of the angular momentum of the rotating object. You can use CEL expressions to calculate the aero/hydrodynamic torque on the object, include any external torques acting on it, and create a Source for the angular momentum additional variable that represents the net torque. Finally you use CEL expressions to physically rotate the mesh.

Note that this only works for transient runs, and that you do NOT use Domain Motion > Rotating and specify an angular velocity; instead, you are physically rotating the mesh.

[ghorrocks]

Not sure about this, but can't you do this without the additional variable and moving mesh? Simply a CEL callback function to read the angular velocity (eg ave(Angular Velocity)@RotatingDomain)? Callback functions are evaluated at the start of a timestep so in effect give you the angular velocity of the last timestep. Then you can get the torques from the model and from external influences and work out the angular velocity of the next timestep. And I understand this works for rotating frame of reference too.

But I admit I have not tried it, just heard of it being done that way. Hopefully it works, it would be much simpler, easier and more efficient than the method you describe.

[michael_owen]

Originally Posted by ghorrocks

Not sure about this, but can't you do this without the additional variable and moving mesh? Simply a CEL callback function to read the angular velocity (eg ave(Angular Velocity)@RotatingDomain)? Callback functions are evaluated at the start of a timestep so in effect give you the angular velocity of the last timestep. Then you can get the torques from the model and from external influences and work out the angular velocity of the next timestep. And I understand this works for rotating frame of reference too.

But I admit I have not tried it, just heard of it being done that way. Hopefully it works, it would be much simpler, easier and more efficient than the method you describe.

Glenn,

I'm afraid that won't work on a couple of accounts. One, the angular velocity of a rotating domain cannot be a function of any variable except time, and you cannot call an integrative function of a non-field variable.

[singer1812]

What Glenn describes is correct. This can be done entirely within CEL as he describes. But there are advantages to putting this in as an additional variable. In the event your motion doesn't simply rotate a subdomain, but instead it moves the mesh, in the event of a remesh incident, you will need the information of angular velocity to be interpolated onto the new mesh file so that you can correct keep the momentum.

[michael_owen]

Try it then.

[singer1812]

Have already. Done it often.

[ghorrocks]

Can you post the CCL of an example?

[singer1812]

Glenn:

Which case are you refering to? CEL or AV case.

[ghorrocks]

The CEL case.

[michael_owen]

I think that there are certain cases where it can be done entirely with CEL without using an additional variable. If you have a surface of constant radius on which you can evalute areaAve(Mesh Velocity), you can turn this into the angular velocity, and then reference that in the CEL expressions for defining the mesh motion.

But again, this involves the method I already described, physically moving the mesh via CEL, and does NOT involve the use of Domain Motion > Rotating, which I repeat will NOT work, because the angular velocity of a rotating frame of reference can only be a function of time.

[ghorrocks]

You can get the angular velocity with something like maxVal(Mesh Velocity)@rotatingdomain where the maximum radius point is known. This can be worked back to the angular velocity for any geometry.

[michael_owen]

Originally Posted by ghorrocks

You can get the angular velocity with something like maxVal(Mesh Velocity)@rotatingdomain where the maximum radius point is known. This can be worked back to the angular velocity for any geometry.

Yeah, I agree. The AV is definitely unnecessary. The mesh motion is not.

[singer1812]

I see where the confusion is. I glossed over the Domain Motion > Rotating statement. I concur with what has been talked about. For thread completeness, here is the CEL you inquired about. This has wall motion rotating around X0,Y0,Z0, along a hinge point pointing in the global Z direction, but not coincedent. The new displacements are updated using New X, New Y and New Z.


LIBRARY:
CEL:
&replace EXPRESSIONS:
Ifin = 127.00 [kg mm^2]
New X = ((x-X0)-Total Mesh Displacement X)*cos(alpha)+((y-Y0)-Total Mesh Displacement Y)*sin(alpha)+X0
New Y = -((x-X0)-Total Mesh Displacement X)*sin(alpha)+((y-Y0)-Total Mesh Displacement Y)*cos(alpha)+Y0
New Z = ((z)-Total Mesh Displacement Z) + Z0
X0 = 44.4355 [mm]
Y0 = 25.181 [mm]
Z0 = move
alpha = ((angdisp -(mit+omgold)*tstep*1[s^-1])+alphstart)*1[deg]
angdisp = initang-newang
back max z = minVal(Global Z Coordinate)@Proj Rear Wall
back min z = maxVal(Global Z Coordinate)@Gun Muzzle Exit
back move = move*max(0, (min(1,(back min z-z)/ (back min z- (back max z) ) ) ))
finforx = (forx1+forx2+forx3+forx4+forx5)
finfory = (fory1+fory2+fory3+fory4+fory5)
fintorq = (finforx*Y0-finfory*X0+(tor1+tor2+tor3+tor4+tor5))
fintorqreal = (tor1+tor2+tor3+tor4+tor5)
forx1 = force_x()@Proj Fin Bottom
forx2 = force_x()@Proj Fin Top
forx3 = force_x()@Proj Fin Hinge
forx4 = force_x()@Proj Fin Side Edges
forx5 = force_x()@Proj Fin Outside Edge
fory1 = force_y()@Proj Fin Bottom
fory2 = force_y()@Proj Fin Top
fory3 = force_y()@Proj Fin Hinge
fory4 = force_y()@Proj Fin Side Edges
fory5 = force_y()@Proj Fin Outside Edge
front max z = maxVal(Global Z Coordinate)@Proj Front Wall
front min z = minVal(Global Z Coordinate)@Proj Front Wall
front move = move* (min(1,(front max z-z)/ (front max z- (front min z + 0.1 [m]) ) ) )
initang = atan((ycod-tmdy)/(xcod-tmdx))/1 [deg]+step(-(xcod-tmdx)/1[m]-0.00001)*180
mit = fintorq/Ifin*tstep*57.296*1[s]
newang = atan(ycod/xcod)/1 [deg]+step(-xcod/1[m]-0.00001)*180
omgold = -(((areaAve(Mesh Velocity X)@Proj Fin Outside Edge)^2+(areaAve(Mesh Velocity Y)@Proj Fin Outside Edge)^2)^0.5/((areaAve(Global X Coordinate)@Proj Fin Outside Edge-X0)^2+(areaAve(Global Y Coordinate)@Proj Fin Outside Edge-Y0)^2)^0.5)*57.296*1[s]
periodmove = front move*step((z-minVal(Global Z Coordinate)@Proj Rear Wall)/1[m])+back move*step((maxVal(Global Z Coordinate)@Gun Muzzle Exit-z)/1[m])
pinit = pstart*zstartramp*(xystartstep)*(1-underfinstartp)+underfinstartp*5e7[Pa]
pstart = 4.5e7[Pa]
tempstart = 1422[K]
tend = 10000[s]
tmdx = areaAve(Total Mesh Displacement X)@Proj Fin Outside Edge
tmdy = areaAve(Total Mesh Displacement Y)@Proj Fin Outside Edge
tor1 = torque_z()@Proj Fin Bottom
tor2 = torque_z()@Proj Fin Top
tor3 = torque_z()@Proj Fin Hinge
tor4 = torque_z()@Proj Fin Side Edges
tor5 = torque_z()@Proj Fin Outside Edge
tstart = 0.0 [s]
tstep = 1e-8 [s]
xcod = areaAve(Global X Coordinate)@Proj Fin Outside Edge-X0
xystartstep = step(0.08 - ((x^2+y^2)/1[m^2])^0.5)
ycod = areaAve(Global Y Coordinate)@Proj Fin Outside Edge-Y0
zstartramp = max(0,1-z/0.04[m])
zvstartramp = max(0,1-z/0.0142[m])
END
END
END

[guoxj]

My problem has been solved.Thanks for all the valuable advices