Compiler Explorer

Source code

module mygrid1
    !!   This module contains...
   use, intrinsic :: iso_fortran_env, only: real64
   implicit none
   private

public :: grid

integer, parameter :: rk = real64

type :: grid
        !! 1D grid
      character(20) :: name = "x [-]"
        !! variable name
      real(rk), allocatable :: edges(:)
        !! vector(0:nc) of cell edges
      real(rk), allocatable :: center(:)
        !! vector(nc) of cell centers, \( x_i \)
      real(rk), allocatable :: width(:)
        !! vector(nc) of cell widths,  \( x_{i+1/2} - x_{i-1/2} \)
      real(rk), pointer :: right(:)
        !! vector(nc) of right cell boundaries, , \( x_{i+1/2} \)
      integer :: ncells
        !! number of cells
   contains
      procedure, pass(self) :: new
      procedure, pass(self) :: left
   end type grid

contains

pure subroutine new(self, xmin, xmax, nc)
      ! Constructor grid1
      class(grid), intent(inout), target :: self
        !! object
      real(rk), intent(in) :: xmin
        !! lower boundary of grid domain
      real(rk), intent(in) :: xmax
        !! upper boundary of grid domain
      integer, intent(in) :: nc
        !! number of grid cells

real(rk) :: xedges(0:nc), rx
      integer :: i

! Compute linear mesh
      rx = (xmax - xmin)/nc
      do concurrent(i=0:nc)
         xedges(i) = xmin + rx*i
      end do

! Map values to grid object
      self%ncells = nc
      self%edges = xedges
      self%right => self%edges(1:nc)
      self%center = (self%left([(i, i=1, nc)]) + self%right)/2
      self%width = self%right - self%left([(i, i=1, nc)])

end subroutine new

elemental real(rk) function left(self,celli)
   class(grid), intent(in) :: self
   integer, intent(in) :: celli
   left = self%edges(celli-1)
   end function left

end module mygrid1

program testgrid

use mygrid1, only: grid
   use, intrinsic :: iso_fortran_env, only: real64
   implicit none
   integer, parameter :: rk = real64

type(grid) :: gx
   real(rk) :: xmin, xmax, tmp
   integer :: nc

xmin = 0._rk
   xmax = 100._rk
   nc = 4
   call gx%new(xmin, xmax, nc)
   tmp = gx%left(1)
   print *, gx%left(1)
   tmp = gx%right(1)
   print *, gx%right

end program