make_icosphere Subroutine

public subroutine make_icosphere(mesh, radius, subdivisions, center)

正二十面体の各面を再帰的に 4 分割し、辺中点を球面へ射影した球面三角形メッシュを生成する。 UV 球と異なり極の特異点がなく、要素面積の最大/最小比は分割数によらず 1.3 程度に収まる。 隣接面の共有辺の中点は同じ和 a+b から作るため、頂点はビット一致し隙間を生じない。

Arguments

Type IntentOptional Attributes Name
type(mesh_type), intent(out) :: mesh

生成した球面三角形メッシュ（要素数 20*4**subdivisions、法線は外向き）。
real(kind=dp), intent(in), optional :: radius

球半径 [m]（省略時 0.5）。
integer(kind=i32), intent(in), optional :: subdivisions

再帰分割回数（省略時 2、0..10）。
real(kind=dp), intent(in), optional :: center(3)

球半径 [m]（省略時 0.5）。

Calls

proc~~make_icosphere~~CallsGraph proc~make_icosphere make_icosphere proc~init_mesh init_mesh proc~make_icosphere->proc~init_mesh proc~update_mesh_geometry update_mesh_geometry proc~init_mesh->proc~update_mesh_geometry proc~build_collision_grid build_collision_grid proc~update_mesh_geometry->proc~build_collision_grid proc~cross~2 cross proc~update_mesh_geometry->proc~cross~2 proc~fill_panel_quadrature fill_panel_quadrature proc~update_mesh_geometry->proc~fill_panel_quadrature proc~init_panel_geometry init_panel_geometry proc~update_mesh_geometry->proc~init_panel_geometry proc~cell_id cell_id proc~build_collision_grid->proc~cell_id proc~coord_to_cell coord_to_cell proc~build_collision_grid->proc~coord_to_cell

Called by

proc~~make_icosphere~~CalledByGraph proc~make_icosphere make_icosphere proc~build_one_template build_one_template proc~build_one_template->proc~make_icosphere proc~build_parts_mesh build_parts_mesh proc~build_parts_mesh->proc~build_one_template proc~build_mesh_from_config build_mesh_from_config proc~build_mesh_from_config->proc~build_parts_mesh proc~load_or_init_run_state load_or_init_run_state proc~load_or_init_run_state->proc~build_mesh_from_config program~main main program~main->proc~load_or_init_run_state

Source Code

  subroutine make_icosphere(mesh, radius, subdivisions, center)
    type(mesh_type), intent(out) :: mesh
    real(dp), intent(in), optional :: radius, center(3)
    integer(i32), intent(in), optional :: subdivisions
    real(dp) :: r, c(3), t, ab(3), bc(3), ca(3)
    integer(i32) :: nsub, level, i, nface, k
    integer(i32), parameter :: faces(3, 20) = reshape([ &
                                                      1, 12, 6, 1, 6, 2, 1, 2, 8, 1, 8, 11, 1, 11, 12, &
                                                      2, 6, 10, 6, 12, 5, 12, 11, 3, 11, 8, 7, 8, 2, 9, &
                                                      4, 10, 5, 4, 5, 3, 4, 3, 7, 4, 7, 9, 4, 9, 10, &
                                                      5, 10, 6, 3, 5, 12, 7, 3, 11, 9, 7, 8, 10, 9, 2], [3, 20])
    real(dp) :: verts(3, 12)
    real(dp), allocatable :: tri(:, :, :), next(:, :, :)

    r = 0.5d0; nsub = 2; c = 0.0d0
    if (present(radius)) r = radius
    if (present(subdivisions)) nsub = subdivisions
    if (present(center)) c = center
    if (nsub < 0 .or. nsub > 10) error stop "subdivisions must be in [0, 10]"

    t = 0.5d0*(1.0d0 + sqrt(5.0d0))
    verts = reshape([ &
                    -1.0d0, t, 0.0d0, 1.0d0, t, 0.0d0, -1.0d0, -t, 0.0d0, 1.0d0, -t, 0.0d0, &
                    0.0d0, -1.0d0, t, 0.0d0, 1.0d0, t, 0.0d0, -1.0d0, -t, 0.0d0, 1.0d0, -t, &
                    t, 0.0d0, -1.0d0, t, 0.0d0, 1.0d0, -t, 0.0d0, -1.0d0, -t, 0.0d0, 1.0d0], [3, 12])
    verts = verts/sqrt(1.0d0 + t*t)

    nface = 20
    allocate (tri(3, 3, nface))
    do i = 1, nface
      do k = 1, 3
        tri(:, k, i) = verts(:, faces(k, i))
      end do
    end do

    do level = 1, nsub
      allocate (next(3, 3, 4*nface))
      do i = 1, nface
        ab = unit_midpoint(tri(:, 1, i), tri(:, 2, i))
        bc = unit_midpoint(tri(:, 2, i), tri(:, 3, i))
        ca = unit_midpoint(tri(:, 3, i), tri(:, 1, i))
        next(:, :, 4*i - 3) = reshape([tri(:, 1, i), ab, ca], [3, 3])
        next(:, :, 4*i - 2) = reshape([ab, tri(:, 2, i), bc], [3, 3])
        next(:, :, 4*i - 1) = reshape([ca, bc, tri(:, 3, i)], [3, 3])
        next(:, :, 4*i) = reshape([ab, bc, ca], [3, 3])
      end do
      call move_alloc(next, tri)
      nface = 4*nface
    end do

    do i = 1, nface
      do k = 1, 3
        tri(:, k, i) = c + r*tri(:, k, i)
      end do
    end do
    call init_mesh(mesh, tri(:, 1, :), tri(:, 2, :), tri(:, 3, :))
  end subroutine make_icosphere