Calculation of Velocity in Post-Processing

edited April 29 in DualSPHysics v4.4
Does anyone by chance have an example considering how to calculate this properly?

I am trying to do this on the example with Periodicity in Main at location x = 0.2, y = 0 and z = 0.002, using this code:

set dirout2=%dirout%\measuretool
%measuretool% -dirin %diroutdata% -points points_cyl.txt -savecsv %dirout2%/VelPointBottom2 -vars:-all,+vel -onlytype:-all,+fluid
if not "%ERRORLEVEL%" == "0" goto fail

The result from the measuretool seems a bit weird when compared with my own code:

My own code is basically done in Julia as:

h = convert(Float32,0.002449) #m, smoothing length for Periodicty, example 2 in Main
x = convert(Float32,0.2) #y = 0 z = 0 in theory now
#\alphaD from Runout file
αD = 92833.429688 #7/(4*pi*(h/2)^2) #Wendland Quintic 2D
#Abs since a particle on both sides of x should be counted in bool later
q(ra::Float32,rb::Array{Float32,1},h::Float32) = abs.(ra.-rb)/h
#Wendland kernel as in Periodicity
W(q) = αD*(1 .- q/2).^4 .* (2 .* q .+ 1)
#Reading points data for particles and velocity
pos = readVtkArray("PartFluid_",PostSPH.Points)
vel = readVtkArray("PartFluid_",PostSPH.Vel)
#Initializing x-velocity array, Vab
Vab = []
for i = 1:length(pos)
qab = q(x,pos[i][:,1],h)
#Bool is to ensure on values between 0 and 2 are taken into account
bool = convert(Array{Int},0 .< qab .< 2)
qabF = q(x,pos[i][:,1],h).*bool
Wab = W(qabF)

Maybe some one can help me out / show me an example of why I am wrong?

Kind regards


  • Has no one else tried to calculate manually and verifying everthing is functioning correctly?

    Kind regards
  • The best way to check it is comparing with reference (analytical, experimental) velocity values of an experiment... as we have done in some papers already

    And compute velocity far from boundary layer since maybe you are not considering GAP fluid-bound or other boundary effects

  • Thanks for your comment. Will try to make the code work for some random point in the domain and see if I can get results, corresponding to DualSPHysics own tool.

    Kind regards
  • Note that we also:
    1) normalised the kernel
    2) we only compute variables if the sum_wab is higher than a minimum value (0.5 for example)... this helps to avoid computing variables when we have only 2-3 particles in the kernel domain
  • that is why we use "kclimit" to define the minimum sum_wab (being maxim 1) to normalise, but also to define a dummy value when sum_wab is not big enough (like less than 0.5 for example)
    dummy values are "0" for velocity, but can be -999999999 for pressure, etc

  • I will try to see if I can find the specific code you used, it must be in one of the files, but thanks so much for your explanation. I am trying to increase the functionality of my post-processing tool,, so it is nice getting some help.

    Will let you know if I couldn't find it, in the weekend.

    Kind regards

    "Kernel correction is also applied when the summation of the kernel values around the position is higher than a value (-kclimit:) defining a dummy value if the correction is not applied (-kcdummy:). "
Sign In or Register to comment.