output_overall_vtu.F90 Source File


This file depends on

sourcefile~~output_overall_vtu.f90~~EfferentGraph sourcefile~output_overall_vtu.f90 output_overall_vtu.F90 sourcefile~output_interface.f90 output_interface.F90 sourcefile~output_overall_vtu.f90->sourcefile~output_interface.f90 sourcefile~control.f90 control.F90 sourcefile~output_interface.f90->sourcefile~control.f90 sourcefile~core.f90 core.F90 sourcefile~output_interface.f90->sourcefile~core.f90 sourcefile~domain.f90 domain.F90 sourcefile~output_interface.f90->sourcefile~domain.f90 sourcefile~input.f90 input.F90 sourcefile~output_interface.f90->sourcefile~input.f90 sourcefile~matrix.f90 matrix.F90 sourcefile~output_interface.f90->sourcefile~matrix.f90 sourcefile~project_settings.f90 project_settings.F90 sourcefile~output_interface.f90->sourcefile~project_settings.f90 sourcefile~properties.f90 properties.F90 sourcefile~output_interface.f90->sourcefile~properties.f90 sourcefile~iteration.f90 iteration.F90 sourcefile~control.f90->sourcefile~iteration.f90 sourcefile~time.f90 time.F90 sourcefile~control.f90->sourcefile~time.f90 sourcefile~allocate.f90 allocate.F90 sourcefile~core.f90->sourcefile~allocate.f90 sourcefile~check_range.f90 check_range.F90 sourcefile~core.f90->sourcefile~check_range.f90 sourcefile~deallocate.f90 deallocate.F90 sourcefile~core.f90->sourcefile~deallocate.f90 sourcefile~error.f90 error.F90 sourcefile~core.f90->sourcefile~error.f90 sourcefile~fortran_utils.f90 fortran_utils.F90 sourcefile~core.f90->sourcefile~fortran_utils.f90 sourcefile~string_utils.f90 string_utils.F90 sourcefile~core.f90->sourcefile~string_utils.f90 sourcefile~types.f90 types.F90 sourcefile~core.f90->sourcefile~types.f90 sourcefile~unique.f90 unique.F90 sourcefile~core.f90->sourcefile~unique.f90 sourcefile~vtk.f90 vtk.F90 sourcefile~core.f90->sourcefile~vtk.f90 sourcefile~vtk_constants.f90 vtk_constants.F90 sourcefile~core.f90->sourcefile~vtk_constants.f90 sourcefile~adjacency.f90 adjacency.F90 sourcefile~domain.f90->sourcefile~adjacency.f90 sourcefile~domain_manager.f90 domain_manager.F90 sourcefile~domain.f90->sourcefile~domain_manager.f90 sourcefile~element.f90 element.F90 sourcefile~domain.f90->sourcefile~element.f90 sourcefile~element_factory.f90 element_factory.F90 sourcefile~domain.f90->sourcefile~element_factory.f90 sourcefile~multicoloring.f90 multicoloring.F90 sourcefile~domain.f90->sourcefile~multicoloring.f90 sourcefile~reordering.f90 reordering.F90 sourcefile~domain.f90->sourcefile~reordering.f90 sourcefile~side.f90 side.F90 sourcefile~domain.f90->sourcefile~side.f90 sourcefile~side_factory.f90 side_factory.F90 sourcefile~domain.f90->sourcefile~side_factory.f90 sourcefile~input_interface.f90 input_interface.F90 sourcefile~input.f90->sourcefile~input_interface.f90 sourcefile~matrix_base.f90 matrix_base.F90 sourcefile~matrix.f90->sourcefile~matrix_base.f90 sourcefile~matrix_coo.f90 matrix_coo.F90 sourcefile~matrix.f90->sourcefile~matrix_coo.f90 sourcefile~matrix_crs.f90 matrix_crs.F90 sourcefile~matrix.f90->sourcefile~matrix_crs.f90 sourcefile~matrix_dense.f90 matrix_dense.F90 sourcefile~matrix.f90->sourcefile~matrix_dense.f90 sourcefile~project_settings.f90->sourcefile~core.f90 sourcefile~materials_manager.f90 materials_manager.F90 sourcefile~properties.f90->sourcefile~materials_manager.f90 sourcefile~properties_manager.f90 properties_manager.F90 sourcefile~properties.f90->sourcefile~properties_manager.f90 sourcefile~adjacency_element.f90 adjacency_element.F90 sourcefile~adjacency.f90->sourcefile~adjacency_element.f90 sourcefile~adjacency_node.f90 adjacency_node.F90 sourcefile~adjacency.f90->sourcefile~adjacency_node.f90 sourcefile~allocate.f90->sourcefile~error.f90 sourcefile~deallocate.f90->sourcefile~error.f90 sourcefile~domain_manager.f90->sourcefile~core.f90 sourcefile~domain_manager.f90->sourcefile~input.f90 sourcefile~domain_manager.f90->sourcefile~adjacency.f90 sourcefile~domain_manager.f90->sourcefile~element.f90 sourcefile~domain_manager.f90->sourcefile~element_factory.f90 sourcefile~domain_manager.f90->sourcefile~multicoloring.f90 sourcefile~domain_manager.f90->sourcefile~reordering.f90 sourcefile~domain_manager.f90->sourcefile~side.f90 sourcefile~domain_manager.f90->sourcefile~side_factory.f90 sourcefile~element.f90->sourcefile~core.f90 sourcefile~element.f90->sourcefile~input.f90 sourcefile~element_factory.f90->sourcefile~core.f90 sourcefile~element_factory.f90->sourcefile~input.f90 sourcefile~element_factory.f90->sourcefile~element.f90 sourcefile~memory_stats_wrapper.f90 memory_stats_wrapper.F90 sourcefile~fortran_utils.f90->sourcefile~memory_stats_wrapper.f90 sourcefile~signal_flag_wrapper.f90 signal_flag_wrapper.F90 sourcefile~fortran_utils.f90->sourcefile~signal_flag_wrapper.f90 sourcefile~system_info_wrapper.f90 system_info_wrapper.F90 sourcefile~fortran_utils.f90->sourcefile~system_info_wrapper.f90 sourcefile~input_interface.f90->sourcefile~core.f90 sourcefile~input_interface.f90->sourcefile~project_settings.f90 sourcefile~materials_manager.f90->sourcefile~input.f90 sourcefile~calculate.f90 calculate.F90 sourcefile~materials_manager.f90->sourcefile~calculate.f90 sourcefile~matrix_base.f90->sourcefile~domain.f90 sourcefile~matrix_coo.f90->sourcefile~core.f90 sourcefile~matrix_coo.f90->sourcefile~domain.f90 sourcefile~matrix_coo.f90->sourcefile~matrix_base.f90 sourcefile~matrix_crs.f90->sourcefile~core.f90 sourcefile~matrix_crs.f90->sourcefile~domain.f90 sourcefile~matrix_crs.f90->sourcefile~matrix_base.f90 sourcefile~matrix_crs.f90->sourcefile~matrix_coo.f90 sourcefile~matrix_dense.f90->sourcefile~core.f90 sourcefile~matrix_dense.f90->sourcefile~domain.f90 sourcefile~matrix_dense.f90->sourcefile~matrix_base.f90 sourcefile~multicoloring.f90->sourcefile~core.f90 sourcefile~multicoloring.f90->sourcefile~adjacency_element.f90 sourcefile~properties_manager.f90->sourcefile~core.f90 sourcefile~properties_manager.f90->sourcefile~input.f90 sourcefile~properties_manager.f90->sourcefile~materials_manager.f90 sourcefile~properties_manager.f90->sourcefile~calculate.f90 sourcefile~reordering.f90->sourcefile~core.f90 sourcefile~reordering.f90->sourcefile~element.f90 sourcefile~reordering.f90->sourcefile~adjacency_node.f90 sourcefile~side.f90->sourcefile~core.f90 sourcefile~side.f90->sourcefile~input.f90 sourcefile~side_factory.f90->sourcefile~core.f90 sourcefile~side_factory.f90->sourcefile~input.f90 sourcefile~side_factory.f90->sourcefile~side.f90 sourcefile~string_utils.f90->sourcefile~allocate.f90 sourcefile~time.f90->sourcefile~core.f90 sourcefile~time.f90->sourcefile~input.f90 sourcefile~array.f90 array.F90 sourcefile~types.f90->sourcefile~array.f90 sourcefile~gauss.f90 gauss.F90 sourcefile~types.f90->sourcefile~gauss.f90 sourcefile~pointer.f90 pointer.F90 sourcefile~types.f90->sourcefile~pointer.f90 sourcefile~variable.f90 variable.F90 sourcefile~types.f90->sourcefile~variable.f90 sourcefile~vector.f90 vector.F90 sourcefile~types.f90->sourcefile~vector.f90 sourcefile~unique.f90->sourcefile~allocate.f90 sourcefile~vtk.f90->sourcefile~allocate.f90 sourcefile~vtk.f90->sourcefile~deallocate.f90 sourcefile~vtk.f90->sourcefile~unique.f90 sourcefile~vtk.f90->sourcefile~vtk_constants.f90 sourcefile~vtk.f90->sourcefile~array.f90 sourcefile~vtk_wrapper.f90 vtk_wrapper.F90 sourcefile~vtk.f90->sourcefile~vtk_wrapper.f90 sourcefile~vtu_wrapper.f90 vtu_wrapper.F90 sourcefile~vtk.f90->sourcefile~vtu_wrapper.f90 sourcefile~adjacency_element.f90->sourcefile~core.f90 sourcefile~adjacency_element.f90->sourcefile~element.f90 sourcefile~adjacency_node.f90->sourcefile~core.f90 sourcefile~array.f90->sourcefile~allocate.f90 sourcefile~array.f90->sourcefile~deallocate.f90 sourcefile~density_interface.f90 density_interface.F90 sourcefile~calculate.f90->sourcefile~density_interface.f90 sourcefile~gcc_interface.f90 gcc_interface.F90 sourcefile~calculate.f90->sourcefile~gcc_interface.f90 sourcefile~heat_capacity_interface.f90 heat_capacity_interface.F90 sourcefile~calculate.f90->sourcefile~heat_capacity_interface.f90 sourcefile~specific_heat_interface.f90 specific_heat_interface.F90 sourcefile~calculate.f90->sourcefile~specific_heat_interface.f90 sourcefile~thermal_conductivity_interface.f90 thermal_conductivity_interface.F90 sourcefile~calculate.f90->sourcefile~thermal_conductivity_interface.f90 sourcefile~c_utils.f90 c_utils.F90 sourcefile~memory_stats_wrapper.f90->sourcefile~c_utils.f90 sourcefile~signal_flag.f90 signal_flag.F90 sourcefile~signal_flag_wrapper.f90->sourcefile~signal_flag.f90 sourcefile~system_info_wrapper.f90->sourcefile~c_utils.f90 sourcefile~variable.f90->sourcefile~allocate.f90 sourcefile~c_utils.f90->sourcefile~signal_flag.f90 sourcefile~memory_stats.f90 memory_stats.F90 sourcefile~c_utils.f90->sourcefile~memory_stats.f90 sourcefile~system_info.f90 system_info.F90 sourcefile~c_utils.f90->sourcefile~system_info.f90 sourcefile~density_interface.f90->sourcefile~core.f90 sourcefile~density_interface.f90->sourcefile~input.f90 sourcefile~gcc_interface.f90->sourcefile~input.f90 sourcefile~heat_capacity_interface.f90->sourcefile~core.f90 sourcefile~heat_capacity_interface.f90->sourcefile~input_interface.f90 sourcefile~heat_capacity_interface.f90->sourcefile~density_interface.f90 sourcefile~specific_heat_interface.f90->sourcefile~core.f90 sourcefile~specific_heat_interface.f90->sourcefile~input.f90 sourcefile~thermal_conductivity_interface.f90->sourcefile~core.f90 sourcefile~thermal_conductivity_interface.f90->sourcefile~input.f90

Source Code

submodule(input_output) input_output_overall_vtu
    implicit none

contains
    module subroutine initialize_output_overall_vtu(self, Input, Coordinate, Domain)
        implicit none
        class(type_output_overall), intent(inout) :: self
        type(Type_Input), intent(in) :: Input
        type(type_dp_3d), intent(in) :: Coordinate
        type(type_domain), intent(inout) :: Domain

        integer(int32) :: i, j
        integer(int32) :: total_connectivity_size, current_offset, start_index

        self%vtk%num_points = input%geometry%vtk%num_points
        self%vtk%num_cells = input%geometry%vtk%num_total_cells
        call self%vtk%coordinate%initialize(self%vtk%num_points)
        self%vtk%coordinate = input%geometry%vtk%POINTS

        ! --- offset と CellType 配列を確保 ---
        call allocate_array(self%vtk%offsets, self%vtk%num_cells)
        call allocate_array(self%vtk%cell_types, self%vtk%num_cells)

        ! ★★★ 修正箇所: offset配列を累積和として正しく計算 ★★★
        current_offset = 0
        do i = 1, self%vtk%num_cells
            self%vtk%cell_types(i) = input%geometry%vtk%CELLS(i)%cell_type
            current_offset = current_offset + input%geometry%vtk%CELLS(i)%num_nodes_in_cell
            self%vtk%offsets(i) = current_offset
        end do

        ! --- connectivity配列を正しい合計サイズで確保 ---
        if (self%vtk%num_cells > 0) then
            total_connectivity_size = self%vtk%offsets(self%vtk%num_cells)
        else
            total_connectivity_size = 0
        end if
        call allocate_array(self%VTK%connectivities, total_connectivity_size)

        ! ★★★ 修正箇所: 正しいオフセットを用いてconnectivity配列にデータを格納 ★★★
        do i = 1, self%vtk%num_cells
            ! 各セルの書き込み開始インデックスを計算
            if (i == 1) then
                start_index = 0
            else
                start_index = self%vtk%offsets(i - 1)
            end if

            ! 各セルの接続情報をコピー (VTKは0-based indexなので-1する)
            do j = 1, input%geometry%vtk%CELLS(i)%num_nodes_in_cell
                self%VTK%connectivities(start_index + j) = input%geometry%vtk%CELLS(i)%connectivity(j) - 1
            end do
        end do

        if (associated(self%write_fields)) nullify (self%write_fields)
        self%write_fields => output_overall_vtu_fields

        if (associated(self%write_cell)) nullify (self%write_cell)
        self%write_cell => output_overall_vtu_cell

    end subroutine initialize_output_overall_vtu

    subroutine output_overall_vtu_fields(self, file_counts, domain, porosity, temperature, si, pressure, water_flux)
        implicit none
        class(type_output_overall), intent(inout) :: self
        integer(int32), intent(in) :: file_counts
        type(type_domain), intent(in) :: domain
        real(real64), intent(in), optional :: porosity(:)
        real(real64), intent(in), optional :: temperature(:)
        real(real64), intent(in), optional :: si(:)
        real(real64), intent(in), optional :: pressure(:)
        type(type_dp_3d), intent(in), optional :: water_flux

        type(vtk_file) :: vtu
        integer(int32) :: status
        integer(int32) :: unit_num

        real(real64), allocatable :: original(:), original_vector(:, :)

        integer(int32) :: nsize, i

        character(256) :: output_name

        ! Initialize VTK file
        write (output_name, self%format_output) trim(self%dir_output_field), "Out_", file_counts, self%file_extension

        status = vtu%initialize(format='ascii', filename=trim(output_name), mesh_topology='UnstructuredGrid')

        ! Write data
        status = vtu%xml_writer%write_piece(np=self%vtk%num_points, &
                                            nc=self%vtk%num_cells)
        status = vtu%xml_writer%write_geo(np=self%vtk%num_points, &
                                          nc=self%vtk%num_cells, &
                                          x=self%vtk%coordinate%x, &
                                          y=self%vtk%coordinate%y, &
                                          z=self%vtk%coordinate%z)
        status = vtu%xml_writer%write_connectivity(nc=self%vtk%num_cells, &
                                                   connectivity=self%VTK%connectivities, &
                                                   offset=self%vtk%offsets, &
                                                   cell_type=self%vtk%cell_types)

        ! --- データセクション ---

        do i = 1, size(self%variable_names)
            if (i == 1) status = vtu%xml_writer%write_dataarray(location='node', action='open')
            select case (self%variable_names(i))
            case ("temperature")
                if (present(temperature)) then
                    call allocate_array(original, self%vtk%num_points)
                    call domain%reordering%to_original_value(temperature, original)
                    status = vtu%xml_writer%write_dataarray(data_name='Temperature', &
                                                            x=original)
                    call deallocate_array(original)
                end if
            case ("ice_saturation")
                if (present(si)) then
                    call allocate_array(original, self%vtk%num_points)
                    call domain%reordering%to_original_value(si, original)
                    status = vtu%xml_writer%write_dataarray(data_name='Si', &
                                                            x=original)
                    call deallocate_array(original)
                end if
            case ("thermal_conductivity")
                print *, "Warning: 'thermal_conductivity' is not implemented in VTK output."
            case ("volumetric_heat_capacity")
                print *, "Warning: 'volumetric_heat_capacity' is not implemented in VTK output."
            case ("pressure")
                if (present(pressure)) then
                    call allocate_array(original, self%vtk%num_points)
                    call domain%reordering%to_original_value(pressure, original)
                    status = vtu%xml_writer%write_dataarray(data_name='Pressure', &
                                                            x=original)
                    call deallocate_array(original)
                end if
            case ("water_flux")
                if (present(water_flux)) then
                    call allocate_array(original_vector, 3_int32, self%vtk%num_points)
                    call domain%reordering%to_original_value(water_flux%x, original_vector(:, 1))
                    call domain%reordering%to_original_value(water_flux%y, original_vector(:, 2))
                    call domain%reordering%to_original_value(water_flux%z, original_vector(:, 3))
                    status = vtu%xml_writer%write_dataarray(data_name='waterFlux', &
                                                            x=original_vector(:, 1), &
                                                            y=original_vector(:, 2), &
                                                            z=original_vector(:, 3))
                    call deallocate_array(original_vector)
                end if
            case ("hydraulic_conductivity")
                print *, "Warning: 'hydraulic_conductivity' is not implemented in VTK output."
            end select

        end do
        status = vtu%xml_writer%write_dataarray(location='node', action='close')
        status = vtu%xml_writer%write_piece()

        ! Finalize VTK file
        status = vtu%finalize()

    end subroutine output_overall_vtu_fields

    subroutine output_overall_vtu_cell(self, file_name, variable_name, variable)
        implicit none
        class(type_output_overall), intent(inout) :: self
        character(*), intent(in) :: file_name
        character(*), intent(in) :: variable_name
        integer(int32), intent(in) :: variable(:)

        type(vtk_file) :: vtu
        integer(int32) :: status
        integer(int32) :: iN, iE, idx, i

        status = vtu%initialize(format='ascii', filename=trim(self%dir_output_field)//trim(file_name)//trim(self%file_extension), &
                                mesh_topology='UnstructuredGrid')

        ! Write data
        status = vtu%xml_writer%write_piece(np=self%vtk%num_points, &
                                            nc=self%vtk%num_cells)
        status = vtu%xml_writer%write_geo(np=self%vtk%num_points, &
                                          nc=self%vtk%num_cells, &
                                          x=self%vtk%coordinate%x, &
                                          y=self%vtk%coordinate%y, &
                                          z=self%vtk%coordinate%z)
        status = vtu%xml_writer%write_connectivity(nc=self%vtk%num_cells, &
                                                   connectivity=self%VTK%connectivities, &
                                                   offset=self%vtk%offsets, &
                                                   cell_type=self%vtk%cell_types)

        status = vtu%xml_writer%write_dataarray(location='cell', action='open')
        status = vtu%xml_writer%write_dataarray(data_name=variable_name, x=variable)
        status = vtu%xml_writer%write_dataarray(location='cell', action='close')
        status = vtu%xml_writer%write_piece()

        status = vtu%finalize()

    end subroutine output_overall_vtu_cell

end submodule input_output_overall_vtu