View Issue Details [ Jump to Notes ]  [ Issue History ] [ Print ]  
ID  Project  Category  View Status  Date Submitted  Last Update  

0002666  OpenFOAM  Patch  public  20170816 15:39  20170925 11:13  
Reporter  tniemi  
Assigned To  will  
Priority  low  Severity  tweak  Reproducibility  sometimes  
Status  resolved  Resolution  fixed  
Product Version  dev  
Target Version  Fixed in Version  dev  
Summary  0002666: Alternative way to calculate momentum transfer from lagrangian parcels  
Description  (I couldn't add files to https://bugs.openfoam.org/view.php?id=2282, so I opened this) The problem: Currently, the momentum transfer in KinematicParcel.C is calculated as dUTrans=dt*(Feff.Sp()*(Ures.average()  Uc_)), where Ures.average() is the mean of the particle velocity during the time step. With Euler integration Ures.average() is just the linear average of the velocities before and after the integration. This is fine if the velocity change has been roughly linear, ie. the integration time step is smaller than the response time of a parcel. On the other hand, if the velocity change is not linear, the error can be very large. In the testcase in https://bugs.openfoam.org/view.php?id=2282, the error in momentum transfer is ~13% with 150 micron particle and ~0.02s step size. If the particle size is changed to 1 micron, the momentum source is 2385.25. Proposed fix: Calculate the momentum of the parcel before and after the integration. The momentum sink to continuous phase is then the difference between these two values when the contribution of possible noncoupled forces is subtracted.  
Tags  No tags attached.  
Attached Files 

Notes  
tniemi (reporter) 20170816 16:12 
Note: using analytical integration will help a lot, it works well up to 1e12 parcel size or so in the test case. However, I guess that Euler is quite often used for U and the behavior may be surprising. 
will (manager) 20170823 09:51 Last edited: 20170823 09:52 
Sorry for the delay here. I think what you are saying makes sense. I don't think the average of the velocity over the integration step multiplied by the force coefficient is a particularly good way of calculating the momentum change. Beyond that, the average velocity over the Euler integration step is not (u_0 + u_1)/2. Euler has the velocity profile with time as asymptotic, so the average (as determined by integration) is a logarithmic function. This clearly removes the advantage of Euler integration as the "cheap" option. Taking the difference is a much better idea. The complication is that the uncoupled forces should not be included, so your solution subtracts them from the final momentum difference. This is probably OK, so long as there aren't any implicit (Sp term) uncoupled forces. It is currently assumed that this is the case. If there are implicit uncoupled forces, then I don't know how to isolate coupled momentum change from uncoupled. Also, it is worth noting that this mechanism is in use for energy, too. I suppose all the same points apply. 
tniemi (reporter) 20170823 14:13 
> If there are implicit uncoupled forces, then I don't know how to isolate coupled momentum change from uncoupled. Yes, I don't know if there is a robust solution for this without resorting to using average velocity. But as you said, currently there are no models using noncoupled Sp and this is also assumed in the current implementation. > Also, it is worth noting that this mechanism is in use for energy, too. I suppose all the same points apply. Yes, probably so. I haven't been that concerned with heat transfer as I have always used analytical for that, which works ok because it gives the correct average. I would guess that analytical is a more common and reasonable choice here, because with a constant htc it gives the exact solution. This makes the possible performance penalty more justified. Also, the problems might not be that bad with heat transfer. With drag, if there is eg. turbulent dispersion, the slip velocity can be temporarily quite large. Even with not so small parcels this can lead to a huge coupling factor/small timescales. The heat transfer coeff tends to react more mildly to slip velocity. 
tniemi (reporter) 20170823 14:48 
Actually, now that I thought about it, it might be possible to calculate the momentum transfer without knowing the velocity. If one also removes the contribution of coupled Su, the remaining momentum change has to be explained by the Sp terms. By using the ratio of uncoupled/coupled Sp one could calculate the coupled portion: coupled mom = Fcp.Sp/(Fcp.Sp+Fncp.Sp)*(mom without Fcp.su) + Fcp.su*dt 
will (manager) 20170825 10:13 
I don't understand where that expression has come from. Could you elaborate? I have an alternative proposal. All the integration routines can be rewritten as a fully explicit contribution multiplied by a factor, F, that depends only on the implicit coefficient. The explicit contribution can be trivially split into it's coupled and noncoupled components. See the uploaded image for the equations. I can't justify the split physically (i.e., why the explicit contribution splits, but F does not), but it's straightforward, conservative, and it tends to the correct form in all the limiting cases I can think of. 
tniemi (reporter) 20170901 12:20 
Sorry for the long delay. I have attached a derivation for my proposal. Case 1 corresponds to original suggestion, case 2 for nonzero noncoupled Sp. You proposal also looks reasonable for splitting the change in velocity to coupled and noncoupled parts, but I'm not sure how to use it to get the momentum transfer. Maybe I'm just thinking this wrong, but let's consider the following: Assume that there is a particle in gas with a drag force and gravity and that the slip velocity is such that drag and gravity are in balance. Then physically delta up should be 0 and the momentum source should be equal to Sp_c*(ugup)*dt (where ug=gas velocity, up=particle velocity). Now in this case ac=Sp_c*ug and bc=Sp_c Using the split formulation and Euler, du_c=(Sp_c*(ugu0))F=(Sp_c(*ugu0))*[dt/(1+Sp_c* dt)] If dt or Sp_c is large, then dt/(1+Sp_c* dt) = 1/Sp_c and du_c=(Sp_c(ugu0))/Sp_c=(ugu0) This tells correctly that given infinite time, the coupled force tries to undo the velocity difference. However, how to get the momentum source from this? 
will (manager) 20170925 11:13 
I've pushed some changes to dev which I think resolves this. Commits a5806207 and fa3ed3b8. The transfer terms are now calculated from the change in particle property, not from evaluating the rate expressions with an average value over the integration step. The integration methods have been rewritten to allow for the split between coupled and noncoupled, and to take advantage of the simplification of not calculating the step average. The Euler scheme is, as far as I can tell, now fully conservative. The difference between UTrans values between Euler and analytical schemes in the #2282 test case is also not as big as before. This has made me confident enough to push. I'm aware that we didn't fully conclude how to handle the coupled/noncoupled split, though. Please reopen the report if you think the calculation still isn't right. 
Issue History  
Date Modified  Username  Field  Change 

20170816 15:39  tniemi  New Issue  
20170816 15:39  tniemi  File Added: momentum.diff  
20170816 16:12  tniemi  Note Added: 0008597  
20170823 09:51  will  Note Added: 0008624  
20170823 09:52  will  Note Edited: 0008624  View Revisions 
20170823 14:13  tniemi  Note Added: 0008625  
20170823 14:48  tniemi  Note Added: 0008626  
20170825 10:13  will  File Added: splitting.png  
20170825 10:13  will  Note Added: 0008633  
20170901 12:17  tniemi  File Added: derivation.pdf  
20170901 12:20  tniemi  Note Added: 0008651  
20170925 11:13  will  Assigned To  => will 
20170925 11:13  will  Status  new => resolved 
20170925 11:13  will  Resolution  open => fixed 
20170925 11:13  will  Fixed in Version  => dev 
20170925 11:13  will  Note Added: 0008798 