I.1.14 Get_Coordinates_Cells_of_Cells_Multi_Mesh Procedure

The main documentation of the Get_Coordinates_Cells_of_Cells_Multi_Mesh Procedure contains additional explanation of this code listing.

  subroutine Get_Coordinates_CoC_MMesh (Coordinates_Cells_of_Cells, Mesh)

    ! Input variable.
  
    type(Multi_Mesh_type), intent(inout) :: Mesh   ! Mesh object.

    ! Input/Output variable.

    ! Coordinates_Cells_of_Cells BNV.
    type(real,3) :: Coordinates_Cells_of_Cells 

    ! Internal variables.

    type(integer,2) :: Flag_Faces_of_Cells ! Mesh boundary flags.

    !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
  
    ! Verify requirements.
  
    VERIFY(Valid_State(Mesh),5)                       ! Mesh is valid.
    ! Coordinates_Cells_of_Cells is valid.
    VERIFY(Valid_State(Coordinates_Cells_of_Cells),5)  
    ! Coordinates_Cells_of_Cells has correct dimensions.
    VERIFY(SIZE(Coordinates_Cells_of_Cells,1) == Mesh%NDimensions,5) 
    VERIFY(SIZE(Coordinates_Cells_of_Cells,2) == Mesh%NCells_PE,5)
    ! Following is not true for polys.
    VERIFY(SIZE(Coordinates_Cells_of_Cells,3) == Mesh%Faces_per_Cell,5)
 
    ! Initializations.

    if (.not.Initialized(Mesh%Coordinates_Cells_DV)) then
      call Initialize (Mesh%Coordinates_Cells_DV, Mesh%Cell_Structure, 2, &
                       dim1=Mesh%NDimensions)
    end if
    if (.not.Initialized(Mesh%Coordinates_Nodes_of_Cells_CA)) then
      call Initialize (Mesh%Coordinates_Nodes_of_Cells_CA, &
                       Mesh%Nodes_of_Cells_Index, 2, &
                       dim1=Mesh%NDimensions)
    end if
    if (.not.Initialized(Mesh%Coordinates_Cells_of_Cells_CA)) then
      call Initialize (Mesh%Coordinates_Cells_of_Cells_CA, &
                       Mesh%Cells_of_Cells_Index, 2, &
                       dim1=Mesh%NDimensions)
    end if

    ! Gather node coordinates to cells.

    Mesh%Coordinates_Nodes_of_Cells_CA = Mesh%Coordinates_Nodes_DV

    ! Calculate cell center coordinates (= average of the node coordinates).

    call Combine_with_Average (Mesh%Coordinates_Cells_DV, &
                               Mesh%Coordinates_Nodes_of_Cells_CA)

    ! Collect coordinates from "other" cells to cells.

    Mesh%Coordinates_Cells_of_Cells_CA = Mesh%Coordinates_Cells_DV
    Coordinates_Cells_of_Cells = Mesh%Coordinates_Cells_of_Cells_CA

    ! Note that the "other" cell coordinates are now correct throughout
    ! the center of the mesh, but incorrect across a boundary face.
    ! 
    ! For nonorthogonal meshes, the "other" cell coordinates across a
    ! boundary face should never be used.
    ! 
    ! For orthogonal meshes, the "other" cell coordinates across a boundary 
    ! face are set to a value which will allow correct modeling of periodic
    ! boundary conditions. The "other" cell index was already set to the
    ! periodic next cell, so the current coordinate values are for cells on 
    ! the far side of the problem. This is fixed by adding/subtracting the 
    ! problem lengths.
    
    if (Mesh%Orthogonality == "Orthogonal") then
      call Initialize (Flag_Faces_of_Cells, NCells_PE(Mesh), &
                       Mesh%Faces_per_Cell)
      call Get_Flag_Faces_of_Cells (Flag_Faces_of_Cells, Mesh)

      ! Fix the edges, since they wrapped around incorrectly.
      ! This is an assumption that is only valid for ortho cells.

      where (Flag_Faces_of_Cells == 1)      ! Left (-x) Face.
        Coordinates_Cells_of_Cells(1,:,:) = &
          Coordinates_Cells_of_Cells(1,:,:) - &
          Mesh%Lengths(1)
      elsewhere (Flag_Faces_of_Cells == 2) ! Right (+x) Face.
        Coordinates_Cells_of_Cells(1,:,:) = &
          Coordinates_Cells_of_Cells(1,:,:) + &
          Mesh%Lengths(1)
      elsewhere (Flag_Faces_of_Cells == 3) ! Front (-y) Face.
        Coordinates_Cells_of_Cells(2,:,:) = &
          Coordinates_Cells_of_Cells(2,:,:) - &
          Mesh%Lengths(2)
      elsewhere (Flag_Faces_of_Cells == 4) ! Back (+y) Face.
        Coordinates_Cells_of_Cells(2,:,:) = &
          Coordinates_Cells_of_Cells(2,:,:) + &
          Mesh%Lengths(2)
      elsewhere (Flag_Faces_of_Cells == 5) ! Bottom (-z) Face.
        Coordinates_Cells_of_Cells(3,:,:) = &
          Coordinates_Cells_of_Cells(3,:,:) - &
          Mesh%Lengths(3)
      elsewhere (Flag_Faces_of_Cells == 6) ! Top (+z) Face.
        Coordinates_Cells_of_Cells(3,:,:) = &
          Coordinates_Cells_of_Cells(3,:,:) + &
          Mesh%Lengths(3)
      end where
      call Finalize (Flag_Faces_of_Cells)
    end if

    ! Verify guarantees.

    VERIFY(Valid_State(Mesh),5)                ! Mesh is valid.
    ! Coordinates_Cells_of_Cells is valid.
    VERIFY(Valid_State(Coordinates_Cells_of_Cells),5) 
  
    return
  end subroutine Get_Coordinates_CoC_MMesh