[Users] Loop seems to miss grid points (Llama)

[Erik Schnetter]

Severin

A general comment: Your variable pio has only 7 digits of accuracy.
You typed more digits, but since this is a single precision constant,
it is rounded to 7 digits. I recommend using "acos(-1.0d0)" instead.

Your loops look correct. If you want to test your code, then you could
initialize a grid function to zero before your loop, and set it to one
or two inside the loop. This will test whether your loop executes for
all points, and which if branch is taken.

There might be another function called after your loop runs. For
example, interpolation from a neighbouring patch might overwrite some
data. You could output a few grid point values after your loop and
compare with what you see in your images to check this.

-erik


On Mon, Dec 2, 2019 at 10:43 AM Severin Frank
<severin.frank at uni-tuebingen.de> wrote:
>
> Dear all,
>
> I'm calculating radial oscillations on a Neutron Star within a
> multipatch (Llama) environment. The coordinate system is Thornburg04.
> For my time evolution of the fluid element displacement equations, I'm
> scheduling routines for my boundary conditions in MoL_PostStep. The
> conditions are with respect to the (changing) surface of the star
> TOV_surface. I want to loop over all grid points to apply the conditions
> on certain grid points (close to the surface) in Fortran.
>
> !$OMP PARALLEL DO private(i,j,k)
>    do k = 1, cctk_lsh(3)
>       do j = 1, cctk_lsh(2)
>          do i = 1, cctk_lsh(1)
>            if (grid_r(i,j,k) >= (TOV_surface-rprec) .AND. grid_r(i,j,k)
> <= (TOV_surface+ rprec)) then
>            drXi(i,j,k) = 0
>            Pidot(i,j,k) = Wrinv(i,j,k) * P(i,j,k) * drrXi(i,j,k) +
> Qr(i,j,k) * Wrinv(i,j,k) * Xi(i,j,k)
>            else if (grid_r(i,j,k) > (TOV_surface + rprec)) then
>            drXi(i,j,k) = 0
>            drrXi(i,j,k) = 0
>            Xi(i,j,k) = 0
>             Pi(i,j,k) = 0
>            Xidot(i,j,k) = 0
>            Pidot(i,j,k) = 0
>            Xeta(i,j,k) = 0
>            end if
>          end do
>       end do
>    end do
> !$OMP END PARALLEL DO
>
> I attached the specific files of my thorn, the parameter file, as well
> as some screenshots of the evolved data.
>
>
> As you can see, the boundary conditions do not affect any part along the
> x axis. Does any one have an idea why my loop seems to "miss" some if-cases?
>
>
> Thanks a lot and best regards,
>
> Severin
>
>
> _______________________________________________
> Users mailing list
> Users at einsteintoolkit.org
> http://lists.einsteintoolkit.org/mailman/listinfo/users



-- 
Erik Schnetter <schnetter at cct.lsu.edu>
http://www.perimeterinstitute.ca/personal/eschnetter/