[Users] Loop seems to miss grid points (Llama)
Severin Frank
severin.frank at uni-tuebingen.de
Thu Nov 28 13:42:26 CST 2019
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
-------------- next part --------------
A non-text attachment was scrubbed...
Name: yz.png
Type: image/png
Size: 473797 bytes
Desc: not available
Url : http://lists.einsteintoolkit.org/pipermail/users/attachments/20191128/131dd1cf/attachment-0003.png
-------------- next part --------------
A non-text attachment was scrubbed...
Name: xz.png
Type: image/png
Size: 470498 bytes
Desc: not available
Url : http://lists.einsteintoolkit.org/pipermail/users/attachments/20191128/131dd1cf/attachment-0004.png
-------------- next part --------------
A non-text attachment was scrubbed...
Name: xy.png
Type: image/png
Size: 439842 bytes
Desc: not available
Url : http://lists.einsteintoolkit.org/pipermail/users/attachments/20191128/131dd1cf/attachment-0005.png
-------------- next part --------------
# Schedule definitions for thorn SWTNS
# Evolved variables
STORAGE: pressuredensity
STORAGE: metricpots
STORAGE: trafo
STORAGE: paramcheck
STORAGE: scalar[2] density[2]
STORAGE: drscalar[2]
STORAGE: drrscalar[2]
STORAGE: tovradius
STORAGE: scalardot densitydot
STORAGE: tovsurface
STORAGE: fluidscalar
#################################################
# Startup and initial
#################################################
SCHEDULE SWTNS_startup AT startup
{
LANG: C
OPTIONS: meta
} "Register banner with Cactus"
SCHEDULE SWTNS_register_MoL IN MoL_Register
{
LANG: C
OPTIONS: meta
} "Register variables with MoL"
SCHEDULE SWTNS_init AT initial AFTER TOV_C_Exact
{
LANG: Fortran
} "Initialise the system"
#SCHEDULE SWTNS_inittest IN CCTK_PREREGRID
#{
# LANG: Fortran
#} "TEST the system"
#SCHEDULE SWTNS_regridtest AT CCTK_POSTREGRID AFTER SpatialCoordinates
#{
# LANG: Fortran
#} "TEST the system"
##################################################
# Calculate the RHS
##################################################
SCHEDULE SWTNS_calc_rhs IN MoL_CalcRHS
{
LANG: Fortran
} "Calculate the RHS"
#SCHEDULE SWTNS_surface IN MoL_CalcRHS AFTER SWTNS_calc_rhs
#{
# LANG: Fortran
#} "Calculate TOV radius"
##################################################
# Apply the boundary conditions
##################################################
SCHEDULE SWTNS_outerboundary IN MoL_PostStep
{
LANG: Fortran
} "Apply outer boundaries"
SCHEDULE SWTNS_RHS_outerboundary IN MoL_RHSBoundaries
{
LANG: Fortran
} "Apply MoL RHS outer boundaries (TOV_surface)"
SCHEDULE SWTNS_boundaries IN MoL_PostStep AFTER SWTNS_outerboundary
{
LANG: Fortran
OPTIONS: level
SYNC: scalar density drscalar drrscalar scalardot densitydot
} "Select the boundary condition"
SCHEDULE GROUP ApplyBCs AS SWTNS_ApplyBCs IN MoL_PostStep AFTER SWTNS_boundaries
{
} "Apply boundary conditions"
##################################################
# Analysis: Frequencies, speed of sound?
##################################################
#if (recalculate_rhs)
#{
# SCHEDULE SWTNS_calc_rhs AT analysis
# {
# LANG: Fortran
# SYNC: scalardot densitydot
# STORAGE: scalardot densitydot
# TRIGGERS: scalardot densitydot
# } "Calculate the RHS"
#}
#todo
-------------- next part --------------
A non-text attachment was scrubbed...
Name: boundary.F90
Type: text/x-fortran
Size: 3031 bytes
Desc: not available
Url : http://lists.einsteintoolkit.org/pipermail/users/attachments/20191128/131dd1cf/attachment-0003.bin
-------------- next part --------------
A non-text attachment was scrubbed...
Name: init.F90
Type: text/x-fortran
Size: 11127 bytes
Desc: not available
Url : http://lists.einsteintoolkit.org/pipermail/users/attachments/20191128/131dd1cf/attachment-0004.bin
-------------- next part --------------
A non-text attachment was scrubbed...
Name: rhs.F90
Type: text/x-fortran
Size: 2735 bytes
Desc: not available
Url : http://lists.einsteintoolkit.org/pipermail/users/attachments/20191128/131dd1cf/attachment-0005.bin
-------------- next part --------------
#===============================================================
# Testing SWTNS on a static TOV star
#
# Severin Frank
#===============================================================
#################### Thorns ####################################
ActiveThorns = "Time MoL"
ActiveThorns = "CoordBase CartGrid3d Boundary StaticConformal"
ActiveThorns = "SymBase ADMBase TmunuBase HydroBase InitBase ADMCoupling ADMMacros"
ActiveThorns = "IOUtil Formaline"
ActiveThorns = "SpaceMask CoordGauge Constants LocalReduce LocalInterp AEILocalInterp LoopControl
SphericalSurface Vectors"
#AHFinderDirect CoordinatesSymmetry CoordinatesSymmetry
ActiveThorns = "Carpet CarpetInterp"
# CarpetRegrid2 CarpetLib CarpetTracker
ActiveThorns = "CarpetIOASCII CarpetIOScalar CarpetIOHDF5 CarpetIOBasic"
ActiveThorns = "NaNChecker CarpetReduce"
ActiveThorns = "CarpetIOASCII"
#ML_ADMConstraints
ActiveThorns = "Coordinates Interpolate2 CarpetInterp2 GlobalDerivative SWTNS"
#ADMDerivatives
ActiveThorns = "TOVSolver"
ActiveThorns = "GenericFD"
#ActiveThorns = "EOS_Omni"
#ActiveThorns = "GRHydro"
#NewRad ?McLachlan?
#NewRad::z_is_radial = yes
#################### Terminate #################################
Cactus::terminate = "time"
Cactus::cctk_final_time = 50
#################### Time Integration ##########################
Time::timestep_method = "given"
Time::timestep = 0.09
MoL::ODE_Method = "Generic"
MoL::Generic_Type = "RK"
MoL::MoL_Intermediate_Steps = 4
MoL::MoL_Num_Scratch_Levels = 3
#MoL::verbose = "extreme"
#################### Grid Setup ################################
Carpet::domain_from_multipatch = yes
CartGrid3D::type = "multipatch"
CartGrid3D::set_coordinate_ranges_on = "all maps"
Coordinates::coordinate_system = "Thornburg04"
#CoordinatesSymmetry::reflection_x = yes
#CoordinatesSymmetry::reflection_y = yes
#CoordinatesSymmetry::reflection_z = yes
#Coordinates::symmetry = "+xyz octant"
#CoordinatesSymmetry::rotating_90 = yes
Coordinates::h_cartesian = 0.1
Coordinates::h_radial = 0.1
SWTNS::rprec = 0.1
TOVSolver::tovrprec = 0.1
SWTNS::SWTNS_amplitude = 0.0031830989
SWTNS::gridboundary = 2 # ueberfluessig
Coordinates::sphere_inner_radius = 1 #
Coordinates::sphere_outer_radius = 8 #7.5284600798562327 #100
Coordinates::n_angular = 9
Driver::ghost_size = 3 #5
Coordinates::outer_boundary_size = 3
Coordinates::patch_boundary_size = 3
Coordinates::additional_overlap_size = 3
Interpolate2::interpolator_order = 4
Interpolate2::continue_if_selftest_fails = no
Interpolate2::verbose = no
Coordinates::store_inverse_jacobian = "yes"
GenericFD::jacobian_group = "Coordinates::jacobian"
GenericFD::jacobian_derivative_group = "Coordinates::jacobian2"
GenericFD::jacobian_identity_map = 0
TmunuBase::stress_energy_storage = yes
TmunuBase::stress_energy_at_RHS = yes
TmunuBase::timelevels = 1
TmunuBase::prolongation_type = none
HydroBase::timelevels = 1
ADMMacros::spatial_order = 4
ADMBase::metric_type = "physical"
#ML_ADMConstraints::timelevels = 3
SpaceMask::use_mask = "yes"
#################### Initial Data ##############################
InitBase::initial_data_setup_method = "init_some_levels"
ADMBase::initial_data = "tov"
ADMBase::initial_lapse = "tov"
ADMBase::initial_shift = "tov"
ADMBase::initial_dtlapse = "zero"
ADMBase::initial_dtshift = "zero"
TOVSolver::TOV_Rho_Central[0] = 1.62e-3 #1e-3 1.62e-3 3.24e-3 1.24e-3
TOVSolver::TOV_Gamma = 2.0
TOVSolver::TOV_K = 100.0
#TOVSolver::TOV_Num_Radial = 10000
SWTNS::outer_bound ="dirichlet"
#################### Finite Differencing ########################
SummationByParts::order = 2
SummationByParts::use_dissipation = no
SummationByParts::onesided_interpatch_boundaries = no
SummationByParts::onesided_outer_boundaries = yes
SummationByParts::sbp_upwind_deriv = yes
SummationByParts::sbp_1st_deriv = yes
SummationByParts::sbp_2nd_deriv = yes
#################### NaN Checker ################################
NaNChecker::check_every = 1
NaNChecker::action_if_found = "just warn" #"terminate", "just warn", "abort"
NaNChecker::check_vars = "HydroBase::rho HydroBase::press SWTNS::gappa
TOVSolver::TOV_PHI TOVSolver::TOV_R TOVSolver::TOV_mr
SWTNS::TOV_surface SWTNS::PHI SWTNS::drPHI SWTNS::LAMBDA
SWTNS::P SWTNS::drP SWTNS::Wrinv SWTNS::Qr
SWTNS::param_dr SWTNS::grid_r
SWTNS::Xi SWTNS::Xidot SWTNS::Pi SWTNS::Pidot
SWTNS::Xeta SWTNS::drXi SWTNS::drrXi
SWTNS::dxdr SWTNS::dydr SWTNS::dzdr"
#SWTNS::dxdrr SWTNS::dydrr SWTNS::dzdrr"
#################### Output #####################################
IOBasic::outInfo_every = 1
IOBasic::outInfo_vars = "
#HydroBase::rho
#HydroBase::press
#SWTNS::gappa
SWTNS::Xi
SWTNS::Xeta
#SWTNS::drXi
SWTNS::drrXi
#SWTNS::param_dr
SWTNS::TOV_surface
#TOVSolver::TOV_mr
#SWTNS::drPHI
#SWTNS::Xidot
#SWTNS::Pi
#SWTNS::Pidot
#SWTNS::PHI
#SWTNS::LAMBDA
#TOVSolver::TOV_PHI
#TOVSolver::TOV_R
"
#CarpetIOASCII::out1d_vars = "SWTNS::Xi"
#CarpetIOASCII::out1d_z = yes
IOASCII::compact_format = yes
# IOASCII::use_grid_coordinates = yes
#IOASCII::separate_grids = no
# IOASCII::one_file_per_group =yes
IOASCII::out0D_every = 1
IOASCII::out0D_point_x = 0
IOASCII::out0D_point_y = 0
IOASCII::out0D_point_z = 2
IOASCII::out0D_vars ="
grid::coordinates{out_every=1000000000}
SWTNS::Xeta
SWTNS::Xi
SWTNS::TOV_surface
#SWTNS::PHI{out_every=100000000000000}
#SWTNS::LAMBDA{out_every=100000000000000}
SWTNS::Xidot
SWTNS::Pi
SWTNS::Pidot
#HydroBase::rho{out_every=1000000000}
#HydroBase::press{out_every=1000000000}
"
IOHDF5::out_every = 8
IOHDF5::out_vars = "
grid::coordinates{out_every=100000000000}
Coordinates::Seve_in_cube{out_every=100000000000000}
HydroBase::rho{out_every=1000000000000000}
HydroBase::press{out_every=100000000000000}
#TOVSolver::TOV_PHI{out_every=100000000000000}
SWTNS::gappa{out_every=100000000000000}
#SWTNS::Pr{out_every=100000000000000}
SWTNS::Qr{out_every=100000000000000}
SWTNS::P{out_every=100000000000000}
SWTNS::drP{out_every=100000000000000}
SWTNS::Wrinv{out_every=100000000000000}
SWTNS::dxdr{out_every=100000000000000}
SWTNS::dydr{out_every=100000000000000}
SWTNS::dzdr{out_every=100000000000000}
#SWTNS::dxdrr{out_every=100000000000000}
#SWTNS::dydrr{out_every=100000000000000}
#SWTNS::dzdrr{out_every=100000000000000}
SWTNS::Xi
SWTNS::Xeta
SWTNS::Xidot
SWTNS::Pi
SWTNS::Pidot
SWTNS::drXi
SWTNS::drrXi
SWTNS::grid_r{out_every=1000000000000000}
SWTNS::param_dr{out_every=100000000000000}
TOVSolver::TOV_mr{out_every=100000000000000}
SWTNS::PHI{out_every=100000000000000}
SWTNS::drPHI{out_every=100000000000000}
SWTNS::LAMBDA{out_every=100000000000000}
"
More information about the Users
mailing list