sandbox/Antoonvh/test_filaments.c

    Test for filaments.h

    We want to test the initialization of a Lamb_Oseen vorticity field using a varying number of filament sections for its discretization.

    Vortex tube: The induvidual segments can be notable. The color represent the error 

    set xr [0.09:11]
set logscale xy
set grid
set size square
set xlabel 'segments per core diameter'
set ylabel 'error'
plot 'out'
    Use about one segment per core radius (script)

    Use about one segment per core radius (script)

    reset
set xr [0.005: 0.015]
set grid
set size square
set xlabel 'Average vorticity'
set ylabel 'z'
plot 'prof12' u 2:1 t '12 segments', 'prof22' u 2:1 t '22 segements', 'prof43' u 2:1 w l lw 2 t '43 segements'
    Average vorticity (script)

    Average vorticity (script)

    The convergence with the number of segments is limited. This is due to the finite length of the vortex tube and the limited grid resolution.

    #include "grid/multigrid3D.h"
#include "filaments.h"
#include "profile6.h"
#include "view.h"
double xp = 0.1, yp = -0.2;

coord line (double t) {
  coord C;
  C.x = xp;
  C.y = yp;
  C.z = t;
  return C;
}

int main () {
  L0 = 10;
  X0 = Y0 = Z0 = -L0/2;
  init_grid (64);
  vector omega[];
  for (int nsegs = 5; nsegs < 200; nsegs *= sqrt(2)) {
    get_vor_vector (omega, line, -2*L0, 2*L0, nsegs, Lamb_Oseen);
    char fname[99];
    sprintf (fname, "prof%d", nsegs);
    profile ({omega.z}, z, fname);
    double err = 0;
    scalar error[];
    foreach() {
      double val = 0;
      foreach_child()
	val += exp(-sq(x - xp) - sq(y - yp))/(8*pi);
      err += dv()*fabs(val - omega.z[]);
      error[] = (val - omega.z[]);
    }
    printf ("%g %g\n", nsegs/(4*L0), err);
    view (theta = pi/2, psi = pi/2);
    isosurface ("omega.z", exp(-1)/pi, color = "error", min = -0.2, max = 0.2);
    save ("omg.mp4", opt = "-r 3");
  }
}