CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > CFX

The problem of inner iteration when calling the subroutine

Register Blogs Members List Search Today's Posts Mark Forums Read

Reply
 
LinkBack Thread Tools Display Modes
Old   October 14, 2010, 21:29
Default The problem of inner iteration when calling the subroutine
  #1
New Member
 
Join Date: Oct 2010
Posts: 6
Rep Power: 6
guoxj is on a distinguished road
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.
guoxj is offline   Reply With Quote

Old   October 14, 2010, 21:37
Default
  #2
Senior Member
 
Michael P. Owen
Join Date: Mar 2009
Posts: 195
Rep Power: 7
michael_owen is on a distinguished road
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 is offline   Reply With Quote

Old   October 14, 2010, 21:44
Default
  #3
Senior Member
 
Michael P. Owen
Join Date: Mar 2009
Posts: 195
Rep Power: 7
michael_owen is on a distinguished road
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.
michael_owen is offline   Reply With Quote

Old   October 14, 2010, 22:42
Default
  #4
New Member
 
Join Date: Oct 2010
Posts: 6
Rep Power: 6
guoxj is on a distinguished road
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 is offline   Reply With Quote

Old   October 14, 2010, 23:26
Default
  #5
New Member
 
Join Date: Oct 2010
Posts: 6
Rep Power: 6
guoxj is on a distinguished road
Here i apology for posting new topic behind Filippo's
guoxj is offline   Reply With Quote

Old   October 15, 2010, 09:42
Default
  #6
Senior Member
 
Michael P. Owen
Join Date: Mar 2009
Posts: 195
Rep Power: 7
michael_owen is on a distinguished road
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.
michael_owen is offline   Reply With Quote

Old   October 16, 2010, 05:34
Default
  #7
Super Moderator
 
Glenn Horrocks
Join Date: Mar 2009
Location: Sydney, Australia
Posts: 9,551
Rep Power: 76
ghorrocks has a spectacular aura aboutghorrocks has a spectacular aura aboutghorrocks has a spectacular aura about
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.
ghorrocks is offline   Reply With Quote

Old   October 16, 2010, 15:56
Default
  #8
Senior Member
 
Michael P. Owen
Join Date: Mar 2009
Posts: 195
Rep Power: 7
michael_owen is on a distinguished road
Quote:
Originally Posted by ghorrocks View Post
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.
michael_owen is offline   Reply With Quote

Old   October 19, 2010, 17:46
Default
  #9
Senior Member
 
Edmund Singer P.E.
Join Date: Aug 2010
Location: Minneapolis, MN
Posts: 429
Rep Power: 10
singer1812 is on a distinguished road
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.
singer1812 is offline   Reply With Quote

Old   October 19, 2010, 18:13
Default
  #10
Senior Member
 
Michael P. Owen
Join Date: Mar 2009
Posts: 195
Rep Power: 7
michael_owen is on a distinguished road
Try it then.
michael_owen is offline   Reply With Quote

Old   October 19, 2010, 18:15
Default
  #11
Senior Member
 
Edmund Singer P.E.
Join Date: Aug 2010
Location: Minneapolis, MN
Posts: 429
Rep Power: 10
singer1812 is on a distinguished road
Have already. Done it often.
singer1812 is offline   Reply With Quote

Old   October 19, 2010, 18:29
Default
  #12
Super Moderator
 
Glenn Horrocks
Join Date: Mar 2009
Location: Sydney, Australia
Posts: 9,551
Rep Power: 76
ghorrocks has a spectacular aura aboutghorrocks has a spectacular aura aboutghorrocks has a spectacular aura about
Can you post the CCL of an example?
ghorrocks is offline   Reply With Quote

Old   October 19, 2010, 18:30
Default
  #13
Senior Member
 
Edmund Singer P.E.
Join Date: Aug 2010
Location: Minneapolis, MN
Posts: 429
Rep Power: 10
singer1812 is on a distinguished road
Glenn:

Which case are you refering to? CEL or AV case.
singer1812 is offline   Reply With Quote

Old   October 19, 2010, 18:35
Default
  #14
Super Moderator
 
Glenn Horrocks
Join Date: Mar 2009
Location: Sydney, Australia
Posts: 9,551
Rep Power: 76
ghorrocks has a spectacular aura aboutghorrocks has a spectacular aura aboutghorrocks has a spectacular aura about
The CEL case.
ghorrocks is offline   Reply With Quote

Old   October 19, 2010, 19:14
Default
  #15
Senior Member
 
Michael P. Owen
Join Date: Mar 2009
Posts: 195
Rep Power: 7
michael_owen is on a distinguished road
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.
michael_owen is offline   Reply With Quote

Old   October 19, 2010, 20:02
Default
  #16
Super Moderator
 
Glenn Horrocks
Join Date: Mar 2009
Location: Sydney, Australia
Posts: 9,551
Rep Power: 76
ghorrocks has a spectacular aura aboutghorrocks has a spectacular aura aboutghorrocks has a spectacular aura about
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.
ghorrocks is offline   Reply With Quote

Old   October 19, 2010, 20:22
Default
  #17
Senior Member
 
Michael P. Owen
Join Date: Mar 2009
Posts: 195
Rep Power: 7
michael_owen is on a distinguished road
Quote:
Originally Posted by ghorrocks View Post
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.
michael_owen is offline   Reply With Quote

Old   October 20, 2010, 10:36
Default
  #18
Senior Member
 
Edmund Singer P.E.
Join Date: Aug 2010
Location: Minneapolis, MN
Posts: 429
Rep Power: 10
singer1812 is on a distinguished road
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
singer1812 is offline   Reply With Quote

Old   October 23, 2010, 20:37
Default
  #19
New Member
 
Join Date: Oct 2010
Posts: 6
Rep Power: 6
guoxj is on a distinguished road
My problem has been solved.Thanks for all the valuable advices
guoxj is offline   Reply With Quote

Reply

Thread Tools
Display Modes

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are On
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
Errorin subroutine appeared when applying cavitation model pitisrisuk CFX 1 July 2, 2012 03:36
OpenFoam-1.5 on Solaris -- compilation problem calling octreeDataPoint(.) constructor cincaipatron OpenFOAM Installation 9 January 11, 2010 06:37
User subroutine problem Yuqing Feng CFX 0 May 5, 2004 23:39
Periodic flow boundary condition problem sudha FLUENT 3 April 28, 2004 08:40
iteration problem sisilxp FLUENT 0 April 23, 2004 23:36


All times are GMT -4. The time now is 14:03.