
[Sponsors] 
May 27, 2010, 05:49 
making waves in OF 1.6

#1 
Member
YS
Join Date: Jan 2010
Posts: 66
Rep Power: 7 
Hello,
I wrote a custom b.c. named "regularWaveMakerFixedValue" which is used to generate regular waves at the inlet. It can be easily extended to generate other kind of waves. Basicall it applies solid wall (U=0) at the air part of the patch and wave velocities computed from wave theory at the water part of the patch, just the same way, at least I believe so, as what GroovyBC does at the boundary. But as can be seen from the attached pictures, due to the strong vortex generated in the air phase near the free surface within the domain, the wave shape has been seriously corrupted, and soon after quarter of wave period, negative value of VOF appears at this location (that is the big visible unphysical ripple right next to the boundary) which leads to the crash of the solver. Could anybody kindly advise what's going wrong here or what should I do to work it around? many thanks! I've attached the source code of the b.c. as well as my initialization & setup files. b.t.w, I also tried to apply inflow&outflow boundary condition at the air part of the patch but got quite similar results. Last edited by Ya_Squall2010; May 27, 2010 at 17:44. 

May 27, 2010, 12:44 

#2 
Senior Member
Niels Gjoel Jacobsen
Join Date: Mar 2009
Location: Rotterdam, The Netherlands
Posts: 1,592
Rep Power: 24 
Hi
Are you also applying the boundary condition on the VOFratio, alpha? If not, doing so might solve your problem. Best regards, Niels 

May 27, 2010, 17:33 

#3 
Member
YS
Join Date: Jan 2010
Posts: 66
Rep Power: 7 
yes, I applied the bc to both alpha1 and U.


May 28, 2010, 12:19 

#4 
Member
YS
Join Date: Jan 2010
Posts: 66
Rep Power: 7 
so far, I've changed different kind of waves and all got similar error. Could any one who's done wave simulations give me some advice on the right chioce of schemes and parameter setup? Many thanks!


May 29, 2010, 14:58 

#5 
Member
YS
Join Date: Jan 2010
Posts: 66
Rep Power: 7 
now I've tried to force the air phase velocity equal to zero at the end of each time step, and also got same problem. A 'warm' start can not solve it as well.


May 31, 2010, 17:20 

#6 
Member
YS
Join Date: Jan 2010
Posts: 66
Rep Power: 7 
attached is the wave variables time histories. 2nd Stokes wave was used. Nothing strange in the profiles.
export.png 

June 1, 2010, 02:34 

#7 
Senior Member
Niels Gjoel Jacobsen
Join Date: Mar 2009
Location: Rotterdam, The Netherlands
Posts: 1,592
Rep Power: 24 
Hi
I do remember similar problems, however I do not remember what solved it in the end. Have you tried ramping the solution over one or several periods by multiplying with a smooth function? The problem seems to occur where the velocity gradient is the largest. Are you also applying the vertical velocity component from the wave theory on the boundary? Best regards, Niels 

June 1, 2010, 06:24 

#8 
Member
YS
Join Date: Jan 2010
Posts: 66
Rep Power: 7 
Thanks for attention. If you are talking about the "warm" start, that is multiply the wave variables by a damping function in the initial couple of waves, yes, I did try, but no use. As can be seen from the snapshot in my latest post, the solver crashes just when the surface elevation is about to change its moving direction, but the velocity components don't. I wonder if this is due to the special shape of the wave shape being chosen. I shall try another waves next. I did apply both horizontal & vertical velocity components on the boundary.
I considered the possibility that the mesh size is not enough to resolve the wave profiles between two wave crests. Please see attached for a fine grid result, which seems even worse and the solver crashes a little bit earlier. export.png 

June 1, 2010, 06:34 

#9  
Member
YS
Join Date: Jan 2010
Posts: 66
Rep Power: 7 
Quote:


June 1, 2010, 18:11 

#10 
Member
YS
Join Date: Jan 2010
Posts: 66
Rep Power: 7 
update:
Today I tried following wave profile, the small hump at the thought was removed. But the code still cannot proceed after the instant shown. Big nagetive alpha1 appears at the same place as before, and as can be seen from the snapshot, abnormal big velocity was resulted. I was wondering if this is due to the inaccurate solving of the velocity field, therefore, I tried to change the default PISO scheme in interFoam to PIMPLE and hoped accurate velocity field can be predicted by introducing more outer iterations. But the experiments told me that it didn't help. I got pretty much the same error and the only 'improvement' was that the crash happened a little bit late than before. That's all. results of pimple+vof 1.jpg 2.jpg 

June 2, 2010, 07:24 

#11 
Senior Member
Niels Gjoel Jacobsen
Join Date: Mar 2009
Location: Rotterdam, The Netherlands
Posts: 1,592
Rep Power: 24 
Here are my fvSchemes and fvSolution. I have made a run using your parameters and it works in my end using my implementation.
One thought: Are you sure, that you compute the wave number correctly? I get 1.358925 m^1. I have noticed that if omega and k does not fulfill the dispersion relation, the solution typically diverges. Best regards, Niels 

June 2, 2010, 07:29 

#12 
Member
YS
Join Date: Jan 2010
Posts: 66
Rep Power: 7 
Thanks，I shall have a look at my wave number and will post mine soon.


June 2, 2010, 09:23 

#13 
Member
YS
Join Date: Jan 2010
Posts: 66
Rep Power: 7 
Hi, I found my wave number is 1.358925, which is the same as yours. I've attached my implementation of wave maker B.C., is it possible for you to compile it and have a go on your computer to see what you get? Alternatively, you could post your files and I will have a try. Many thanks!
just uncompress the files to .../OpenFOAM1.6/src/finiteVolume/fields/fvPatchFields/derived/, and change the make file accordingly. regularWaveMakerFixedValue.tar.gz 

June 3, 2010, 08:10 

#14 
Senior Member
Niels Gjoel Jacobsen
Join Date: Mar 2009
Location: Rotterdam, The Netherlands
Posts: 1,592
Rep Power: 24 
The only reason I can come up with for your problems is that you for some reason or the other have a force imbalance, which is tried to be corrected through the presence of large velocities, i.e. large shear stresses.
Hence your gamma distribution or momentum flux at the boundary seems not to model what you are looking at. This could be either due to missing pressure gradient at the boundary or the way you determine wet/dry cells. In my implementation I consider "continuous" variation in alpha instead of binary variation. Adding a pressure gradient seems to have limited impact. Good luck, Niels 

September 28, 2010, 05:59 
WaveBC

#15 
New Member
Pal Schmitt
Join Date: Aug 2010
Location: Belfast
Posts: 19
Rep Power: 6 
Hi,
Hope I am not too late here to get some help... I have implemented a waveBC similar to the one shown above, and seem to face the same problem. What puzzles me most is that the same case and BC implemented in groovyBC runs fine. In the picture You can see that the same velocities and alpha conditions are applied, but inside the domain the surface elevation develops different. Ramping etc. has no effect either. Is there something special in groovyBC which is not obvious to understand when setting time dependant values on the Patch? Are waveBC so sensitive that a slightly different value for solving the same equation in groovyBC or in Foam might give wrong results in one and work perfectly in another? Any other idea? Thanks in advance for any hint... Pal 

September 28, 2010, 07:33 

#16 
Senior Member
Niels Gjoel Jacobsen
Join Date: Mar 2009
Location: Rotterdam, The Netherlands
Posts: 1,592
Rep Power: 24 
Hi Pal
You need to give more details. The waveBC implementation linked to above is for infinite deep water, so that might explain your problems, if you are comparing with a groovyBC version on finite water depth. Best regards, Niels 

September 28, 2010, 09:46 
Clearer

#17 
New Member
Pal Schmitt
Join Date: Aug 2010
Location: Belfast
Posts: 19
Rep Power: 6 
Hi,
OK, I will try to be more precise. I implemented and used a simple waveBC in groovyBC, it works fine. alpha1 Code:
type groovyBC; valueExpression "(pos().z<=Depth+A*sin(w*time())) ? 1 : 0"; variables "Depth=1.2;l=1;A=0.05;Endtime=2;rampfactor= ((time()/Endtime)< 1 ) ? (time()/Endtime) : 1;g=vector(0,0,9.81);k=2*pi/l;w=sqrt(k*mag(g));"; Code:
type groovyBC; valueExpression "(pos().z<=Depth+A*sin(w*time())) ? vector( A*w*exp(k*(pos().zDepth))*sin(w*time()), 0,A*w*exp(k*(pos().zDepth))*cos(w*time())) : vector(0,0,0)"; variables "Depth=1.2;l=1;A=0.05;Endtime=2;rampfactor=((time()/Endtime)< 1 ) ? (time()/Endtime) : 1;g=vector(0,0,9.81);k=2*pi/l;w=sqrt(k*mag(g));"; alpha Code:
forAll(result,i) { if (ZCo[i] <= waterLevel_+*amplitude_*sin(omegaTime)) { result[i] = scalar(1); } else { result[i] = scalar(0); } Code:
forAll(result,i) { if ( ZCo[i] <= waterLevel_+*amplitude_*sin(omegaTime) ) { Depthexpress = exp( wavenumber*(ZCo[i]waterLevel_) ); result[i] = vector( amplitude_*omegawave*Depthexpress*sin(omegaTime), 0, amplitude_*omegawave*Depthexpress*cos(omegaTime) ); } else { result[i] = vector(0,0,0); } If I visualize only the patch, the values shown in paraview are the ones created by the BC, or? So, on the one hand the BC seem to be equal, but on the other hand the results are not. Thus the question if there is some trick in groovyBC... Thanks for Your input, Pal Last edited by wavemaster; September 28, 2010 at 11:37. 

September 28, 2010, 11:27 

#18 
Senior Member
Niels Gjoel Jacobsen
Join Date: Mar 2009
Location: Rotterdam, The Netherlands
Posts: 1,592
Rep Power: 24 
The velocity component in the xdirection for the 2 implementations differ by a "".
Best regards, Niels 

September 28, 2010, 11:44 

#19 
New Member
Pal Schmitt
Join Date: Aug 2010
Location: Belfast
Posts: 19
Rep Power: 6 
Sorry, this was a copy error. It did not affect my simulations.
I had just tried waves running in different directions, with groovyBC all is fine and again my brilliant piece of code fails. If I set the amplitude to zero, there is no velocity or alpha movement visible on the patch , but I still get crazy velocities with my stuff... Thanks for looking at it... Pal 

Thread Tools  
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
OF 1.6  Ubuntu 9.10 (64bit)  GLIBCXX_3.4.11 not found  piprus  OpenFOAM Installation  22  February 25, 2010 14:43 
OpenFOAM on MinGW crosscompiler hosted on Linux  allenzhao  OpenFOAM Installation  127  January 30, 2009 20:08 
DxFoam reader update  hjasak  OpenFOAM PostProcessing  69  April 24, 2008 01:24 
Errors running allwmake in OpenFOAM141dev with WM_COMPILE_OPTION%3ddebug  unoder  OpenFOAM Installation  11  January 30, 2008 21:30 
DecomposePar links against liblamso0 with OpenMPI  jens_klostermann  OpenFOAM Bugs  11  June 28, 2007 17:51 