initialize_simple Subroutine

private subroutine initialize_simple(self, num_nodes, elements)

Arguments

Type IntentOptional Attributes Name
class(type_map_node_to_element), intent(inout) :: self
integer(kind=int32), intent(in) :: num_nodes
type(holder_elements), intent(in) :: elements(:)

Calls

proc~~initialize_simple~~CallsGraph proc~initialize_simple initialize_simple proc~get_connectivity abst_mesh%get_connectivity proc~initialize_simple->proc~get_connectivity proc~get_num_nodes~3 abst_mesh%get_num_nodes proc~initialize_simple->proc~get_num_nodes~3

Called by

proc~~initialize_simple~~CalledByGraph proc~initialize_simple initialize_simple proc~initialize_map type_map_node_to_element%initialize_map proc~initialize_map->proc~initialize_simple

Source Code

    subroutine initialize_simple(self, num_nodes, elements)
        class(type_map_node_to_element), intent(inout) :: self
        integer(int32), intent(in) :: num_nodes
        type(holder_elements), intent(in) :: elements(:)

        integer(int32) :: ielem, inode_local, node_id, num_elements, current_size
        integer(int32), allocatable :: temp_list(:)
        integer(int32), dimension(:), pointer :: connectivity

        num_elements = size(elements)
        allocate (self%map_data(num_nodes))

        do ielem = 1, num_elements
            connectivity => elements(ielem)%e%get_connectivity()
            do inode_local = 1, elements(ielem)%e%get_num_nodes()
                node_id = connectivity(inode_local)
                if (node_id > 0 .and. node_id <= num_nodes) then
                    if (allocated(self%map_data(node_id)%ids)) then
                        current_size = size(self%map_data(node_id)%ids)
                        allocate (temp_list(current_size))
                        temp_list = self%map_data(node_id)%ids
                        deallocate (self%map_data(node_id)%ids)
                        allocate (self%map_data(node_id)%ids(current_size + 1))
                        self%map_data(node_id)%ids(1:current_size) = temp_list
                        self%map_data(node_id)%ids(current_size + 1) = ielem
                        deallocate (temp_list)
                    else
                        allocate (self%map_data(node_id)%ids(1))
                        self%map_data(node_id)%ids(1) = ielem
                    end if
                end if
            end do
        end do
    end subroutine initialize_simple