find_first_hit_base_linear Subroutine

public subroutine find_first_hit_base_linear(mesh, p0, p1, seg_min, seg_max, hit, best_t, use_box_filter, box_min, box_max, box_tol, require_elem_inside)

旧実装と同じ線形探索で最初の命中要素を探索する。

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)
real(kind=dp), intent(in) :: seg_min(3)
real(kind=dp), intent(in) :: seg_max(3)
type(hit_info), intent(inout) :: hit
real(kind=dp), intent(inout) :: best_t
logical, intent(in) :: use_box_filter
real(kind=dp), intent(in) :: box_min(3)
real(kind=dp), intent(in) :: box_max(3)
real(kind=dp), intent(in) :: box_tol
logical, intent(in) :: require_elem_inside

Calls

proc~~find_first_hit_base_linear~~CallsGraph proc~find_first_hit_base_linear find_first_hit_base_linear proc~bbox_inside_box bbox_inside_box proc~find_first_hit_base_linear->proc~bbox_inside_box proc~point_inside_box point_inside_box proc~find_first_hit_base_linear->proc~point_inside_box proc~segment_bbox_overlap_precomputed segment_bbox_overlap_precomputed proc~find_first_hit_base_linear->proc~segment_bbox_overlap_precomputed proc~segment_triangle_intersect segment_triangle_intersect 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_base_linear~~CalledByGraph proc~find_first_hit_base_linear find_first_hit_base_linear proc~find_first_hit_base find_first_hit_base proc~find_first_hit_base->proc~find_first_hit_base_linear proc~find_first_hit find_first_hit proc~find_first_hit->proc~find_first_hit_base proc~find_first_hit_periodic2 find_first_hit_periodic2 proc~find_first_hit->proc~find_first_hit_periodic2 proc~find_first_hit_periodic2->proc~find_first_hit_base 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_base_linear( &
    mesh, p0, p1, seg_min, seg_max, hit, best_t, use_box_filter, box_min, box_max, box_tol, require_elem_inside &
    )
    type(mesh_type), intent(in) :: mesh
    real(dp), intent(in) :: p0(3), p1(3), seg_min(3), seg_max(3)
    type(hit_info), intent(inout) :: hit
    real(dp), intent(inout) :: best_t
    logical, intent(in) :: use_box_filter
    real(dp), intent(in) :: box_min(3), box_max(3), box_tol
    logical, intent(in) :: require_elem_inside

    integer(i32) :: i
    logical :: ok
    real(dp) :: t, h(3)

    do i = 1, mesh%nelem
      if (use_box_filter) then
        if (.not. segment_bbox_overlap_precomputed(mesh%bb_min(:, i), mesh%bb_max(:, i), box_min, box_max)) cycle
        if (require_elem_inside) then
          if (.not. bbox_inside_box(mesh%bb_min(:, i), mesh%bb_max(:, i), box_min, box_max, box_tol)) cycle
        end if
      end if
      if (.not. segment_bbox_overlap_precomputed(seg_min, seg_max, mesh%bb_min(:, i), mesh%bb_max(:, i))) cycle
      call segment_triangle_intersect(p0, p1, mesh%v0(:, i), mesh%v1(:, i), mesh%v2(:, i), ok, t, h)
      if (.not. ok) cycle
      if (use_box_filter) then
        if (.not. point_inside_box(h, box_min, box_max, box_tol)) cycle
      end if
      if (t < best_t) then
        best_t = t
        hit%has_hit = .true.
        hit%elem_idx = i
        hit%t = t
        hit%pos = h
        hit%pos_wrapped = h
        hit%image_shift = 0_i32
      end if
    end do
  end subroutine find_first_hit_base_linear