find_first_hit_periodic2 Subroutine

public subroutine find_first_hit_periodic2(mesh, p0, p1, hit, sim, box_min, box_max, require_elem_inside)

periodic2 用に image shift を列挙し、base collision の結果を統合する。

Arguments

Type IntentOptional Attributes Name
type(mesh_type), intent(in) :: mesh
real(kind=dp), intent(in) :: p0(3)
real(kind=dp), intent(in) :: p1(3)
type(hit_info), intent(out) :: hit
type(sim_config), intent(in) :: sim
real(kind=dp), intent(in), optional :: box_min(3)
real(kind=dp), intent(in), optional :: box_max(3)
logical, intent(in), optional :: require_elem_inside

Calls

proc~~find_first_hit_periodic2~~CallsGraph proc~find_first_hit_periodic2 find_first_hit_periodic2 proc~compute_periodic_shift_bounds compute_periodic_shift_bounds proc~find_first_hit_periodic2->proc~compute_periodic_shift_bounds proc~find_first_hit_base find_first_hit_base proc~find_first_hit_periodic2->proc~find_first_hit_base proc~initialize_hit initialize_hit proc~find_first_hit_periodic2->proc~initialize_hit proc~point_inside_box_periodic2 point_inside_box_periodic2 proc~find_first_hit_periodic2->proc~point_inside_box_periodic2 proc~prefer_periodic_candidate prefer_periodic_candidate proc~find_first_hit_periodic2->proc~prefer_periodic_candidate proc~resolve_box_filter_args resolve_box_filter_args proc~find_first_hit_periodic2->proc~resolve_box_filter_args proc~resolve_periodic2_collision_config resolve_periodic2_collision_config proc~find_first_hit_periodic2->proc~resolve_periodic2_collision_config proc~wrap_periodic2_point wrap_periodic2_point proc~find_first_hit_periodic2->proc~wrap_periodic2_point proc~find_first_hit_base->proc~initialize_hit proc~find_first_hit_base->proc~resolve_box_filter_args proc~find_first_hit_base_grid find_first_hit_base_grid proc~find_first_hit_base->proc~find_first_hit_base_grid proc~find_first_hit_base_linear find_first_hit_base_linear proc~find_first_hit_base->proc~find_first_hit_base_linear proc~lower_ascii lower_ascii proc~resolve_periodic2_collision_config->proc~lower_ascii proc~bbox_inside_box bbox_inside_box proc~find_first_hit_base_grid->proc~bbox_inside_box proc~cell_id cell_id proc~find_first_hit_base_grid->proc~cell_id proc~coord_to_cell coord_to_cell proc~find_first_hit_base_grid->proc~coord_to_cell proc~point_inside_box point_inside_box proc~find_first_hit_base_grid->proc~point_inside_box proc~segment_aabb_intersection_t segment_aabb_intersection_t proc~find_first_hit_base_grid->proc~segment_aabb_intersection_t proc~segment_bbox_overlap_precomputed segment_bbox_overlap_precomputed proc~find_first_hit_base_grid->proc~segment_bbox_overlap_precomputed proc~segment_triangle_intersect segment_triangle_intersect proc~find_first_hit_base_grid->proc~segment_triangle_intersect proc~find_first_hit_base_linear->proc~bbox_inside_box proc~find_first_hit_base_linear->proc~point_inside_box proc~find_first_hit_base_linear->proc~segment_bbox_overlap_precomputed proc~find_first_hit_base_linear->proc~segment_triangle_intersect proc~cross cross proc~segment_triangle_intersect->proc~cross

Called by

proc~~find_first_hit_periodic2~~CalledByGraph proc~find_first_hit_periodic2 find_first_hit_periodic2 proc~find_first_hit find_first_hit proc~find_first_hit->proc~find_first_hit_periodic2 proc~sample_photo_raycast_particles sample_photo_raycast_particles proc~sample_photo_raycast_particles->proc~find_first_hit proc~sample_photo_species_state sample_photo_species_state proc~sample_photo_species_state->proc~sample_photo_raycast_particles proc~init_particle_batch_from_config init_particle_batch_from_config proc~init_particle_batch_from_config->proc~sample_photo_species_state

Source Code

  subroutine find_first_hit_periodic2(mesh, p0, p1, hit, sim, box_min, box_max, require_elem_inside)
    type(mesh_type), intent(in) :: mesh
    real(dp), intent(in) :: p0(3), p1(3)
    type(hit_info), intent(out) :: hit
    type(sim_config), intent(in) :: sim
    real(dp), intent(in), optional :: box_min(3), box_max(3)
    logical, intent(in), optional :: require_elem_inside

    integer(i32) :: periodic_axes(2), image_shift(2), nmin(2), nmax(2), iaxis, n1, n2
    real(dp) :: periodic_len(2), shift_vec(3), candidate_pos(3), candidate_wrapped(3)
    real(dp) :: box_min_local(3), box_max_local(3), box_tol
    real(dp) :: shifted_p0(3), shifted_p1(3)
    type(hit_info) :: candidate
    logical :: use_periodic2, use_box_filter, require_inside_elem

    call initialize_hit(hit)
    if (.not. mesh%periodic2_collision_ready) then
      error stop 'periodic2 collision requires prepare_periodic2_collision_mesh before ray queries.'
    end if

    call resolve_periodic2_collision_config(sim, use_periodic2, periodic_axes, periodic_len)
    if (.not. use_periodic2) then
      call find_first_hit_base(mesh, p0, p1, hit, box_min, box_max, require_elem_inside)
      return
    end if

    call resolve_box_filter_args( &
      box_min, box_max, require_elem_inside, use_box_filter, box_min_local, box_max_local, box_tol, require_inside_elem &
      )

    do iaxis = 1, 2
      call compute_periodic_shift_bounds(mesh, p0, p1, periodic_axes(iaxis), periodic_len(iaxis), nmin(iaxis), nmax(iaxis))
      if (nmin(iaxis) > nmax(iaxis)) return
    end do

    do n1 = nmin(1), nmax(1)
      do n2 = nmin(2), nmax(2)
        image_shift = [n1, n2]
        shift_vec = 0.0d0
        shift_vec(periodic_axes(1)) = real(image_shift(1), dp)*periodic_len(1)
        shift_vec(periodic_axes(2)) = real(image_shift(2), dp)*periodic_len(2)
        shifted_p0 = p0 - shift_vec
        shifted_p1 = p1 - shift_vec

        call find_first_hit_base(mesh, shifted_p0, shifted_p1, candidate)
        if (.not. candidate%has_hit) cycle

        candidate_pos = candidate%pos + shift_vec
        candidate_wrapped = candidate_pos
        call wrap_periodic2_point(candidate_wrapped, sim%box_min, periodic_axes, periodic_len)
        if (use_box_filter) then
          if (.not. point_inside_box_periodic2( &
              candidate_wrapped, box_min_local, box_max_local, box_tol, periodic_axes, require_inside_elem)) cycle
        end if

        if (prefer_periodic_candidate(candidate%t, candidate%elem_idx, image_shift, hit)) then
          hit%has_hit = .true.
          hit%elem_idx = candidate%elem_idx
          hit%t = candidate%t
          hit%pos = candidate_pos
          hit%pos_wrapped = candidate_wrapped
          hit%image_shift = image_shift
        end if
      end do
    end do
  end subroutine find_first_hit_periodic2