find_first_hit_base Subroutine

public subroutine find_first_hit_base(mesh, p0, p1, hit, box_min, box_max, 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)
type(hit_info), intent(out) :: hit
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_base~~CallsGraph proc~find_first_hit_base find_first_hit_base 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~initialize_hit initialize_hit proc~find_first_hit_base->proc~initialize_hit proc~resolve_box_filter_args resolve_box_filter_args proc~find_first_hit_base->proc~resolve_box_filter_args 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_base~~CalledByGraph proc~find_first_hit_base find_first_hit_base 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(mesh, p0, p1, hit, 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
    real(dp), intent(in), optional :: box_min(3), box_max(3)
    logical, intent(in), optional :: require_elem_inside

    real(dp) :: d(3), seg_min(3), seg_max(3), best_t
    real(dp) :: box_min_local(3), box_max_local(3), box_tol
    logical :: use_box_filter, require_inside_elem

    call initialize_hit(hit)

    d = p1 - p0
    seg_min = min(p0, p1)
    seg_max = max(p0, p1)
    best_t = huge(1.0d0)

    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 &
      )

    if (mesh%use_collision_grid) then
      call find_first_hit_base_grid( &
        mesh, p0, p1, d, seg_min, seg_max, hit, best_t, &
        use_box_filter, box_min_local, box_max_local, box_tol, require_inside_elem &
        )
    else
      call find_first_hit_base_linear( &
        mesh, p0, p1, seg_min, seg_max, hit, best_t, &
        use_box_filter, box_min_local, box_max_local, box_tol, require_inside_elem &
        )
    end if
  end subroutine find_first_hit_base