VERY GOOD: Imported grid from FortranOpenAcc!!!

A step towards unifying these two; introduced Grid_Mod
from FortranOpenAcc project to here.

On branch master

modified:   Demo_Mod/Fill_In.f90
modified:   Demo_Mod/Main.f90
modified:   Discretize_Mod/On_Sparse_Matrix.f90
modified:   Grid_Mod.f90
new file:   Grid_Mod/Cell_Number.f90
new file:   Grid_Mod/Cells_I_J_K.f90
new file:   Grid_Mod/Create_Grid.f90
new file:   Grid_Mod/Save_Vtk_Debug.f90
modified:   Solvers_Mod/Dense/Cholesky.f90
modified:   Solvers_Mod/Dense/Gauss.f90
modified:   Solvers_Mod/Dense/Ldlt.f90
modified:   Solvers_Mod/Dense/Lu.f90
modified:   Solvers_Mod/Incomplete/Cholesky.f90
modified:   Solvers_Mod/Incomplete/Ldlt.f90
modified:   Solvers_Mod/Incomplete/Lu.f90
modified:   Solvers_Mod/Incomplete/Tflows_Ldlt.f90
modified:   Solvers_Mod/Sparse/Cg_Cholesky_Prec.f90
modified:   Solvers_Mod/Sparse/Cg_Diag_Prec.f90
modified:   Solvers_Mod/Sparse/Cg_No_Prec.f90
modified:   Solvers_Mod/Sparse/Cg_Tflows_Prec.f90
modified:   Sparse_Mod/Create.f90
modified:   makefile_explicit_dependencies

Demo_Mod/Fill_In.f90

@@ -1,16 +1,16 @@
 !==============================================================================!
-  subroutine Demo_Mod_Fill_In(fill_in, grid)
+  subroutine Demo_Mod_Fill_In(fill_in, Grid)
 !------------------------------------------------------------------------------!
   implicit none
 !---------------------------------[Arguments]----------------------------------!
   integer         :: fill_in
-  type(Grid_Type) :: grid
+  type(Grid_Type) :: Grid
 !-----------------------------------[Locals]-----------------------------------!
   type(Sparse_Type) :: A  ! original matrix
   type(Sparse_Type) :: C  ! preconditioning matrix
 !==============================================================================!
 
-  call A % Sparse_Create(grid)
+  call A % Sparse_Create(Grid, singular=.false.)
 
   call IO % Plot_Sparse ("compressed_a",  A)
   call IO % Print_Sparse("Compressed A:", A)

Demo_Mod/Main.f90

@@ -3,14 +3,16 @@ subroutine Demo_Mod_Main
 !------------------------------------------------------------------------------!
   implicit none
 !-----------------------------------[Locals]-----------------------------------!
-  type(Grid_Type)   :: grid
+  type(Grid_Type)   :: Grid
   type(Dense_Type)  :: Ad
   type(Sparse_Type) :: As
   real, allocatable :: x(:)
   real, allocatable :: b(:)
   character(80)     :: dummy
   character(32)     :: arg = ''              ! command line argument
   integer           :: choice, f_in, n_iter, file_unit, test
+  integer           :: nx, ny, nz
+  real              :: lx, ly, lz
   real              :: res
   logical           :: file_exists
 !==============================================================================!
@@ -27,12 +29,12 @@ subroutine Demo_Mod_Main
   if(.not. file_exists) then
 
     print *, "# File 'init.dat' doesn't exist, setting the default values"
-    grid % lx =  1.0
-    grid % ly =  1.0
-    grid % lz =  1.0
-    grid % nx = 10
-    grid % ny = 10
-    grid % nz = 10
+    lx =  1.0
+    ly =  1.0
+    lz =  1.0
+    nx = 10
+    ny = 10
+    nz = 10
 
     f_in   =  1
     n_iter = 10
@@ -46,18 +48,23 @@ subroutine Demo_Mod_Main
          file    = 'init.dat',  &
          action  = 'read')
 
-    read(file_unit, *) dummy, grid % lx
-    read(file_unit, *) dummy, grid % ly
-    read(file_unit, *) dummy, grid % lz
-    read(file_unit, *) dummy, grid % nx
-    read(file_unit, *) dummy, grid % ny
-    read(file_unit, *) dummy, grid % nz
+    read(file_unit, *) dummy, lx
+    read(file_unit, *) dummy, ly
+    read(file_unit, *) dummy, lz
+    read(file_unit, *) dummy, nx
+    read(file_unit, *) dummy, ny
+    read(file_unit, *) dummy, nz
 
     read(file_unit, *) dummy, f_in
     read(file_unit, *) dummy, n_iter
     read(file_unit, *) dummy, res
   end if
 
+  !-----------------------------------!
+  !   Create the computational grid   !
+  !-----------------------------------!
+  call Grid % Create_Grid(lx, ly, lz, nx, ny, nz)
+
   !------------------------------------------------------!
   !   Check if command line argument has been supplied   !
   !------------------------------------------------------!
@@ -106,7 +113,7 @@ subroutine Demo_Mod_Main
   call Foul % Formatted_Write(' # ', 'default',  &
                          'Section 5 - Various Settings', 'bright cyan');
   write(*,'(a46,3i4)') '# 51 - Change grid resolution, currently at: ',  &
-                       grid % nx, grid % ny, grid % nz
+                       Grid % nx, Grid % ny, Grid % nz
   write(*,'(a46,1i4)') '# 52 - Change fill-in level, currently at:   ',  &
                        f_in
   write(*,'(a46,1i4)') '# 53 - Change num iterations, currently at:  ',  &
@@ -141,46 +148,46 @@ subroutine Demo_Mod_Main
       return
 
     case(11)
-        call Solvers_Mod_Gauss(grid, Ad, x, b)
+        call Solvers_Mod_Gauss(Grid, Ad, x, b)
     case(12)
-        call Solvers_Mod_Cholesky(grid, Ad, x, b)
+        call Solvers_Mod_Cholesky(Grid, Ad, x, b)
     case(13)
-        call Solvers_Mod_Ldlt(grid, Ad, x, b)
+        call Solvers_Mod_Ldlt(Grid, Ad, x, b)
     case(14)
-        call Solvers_Mod_Lu(grid, Ad, x, b, GAUSS)
+        call Solvers_Mod_Lu(Grid, Ad, x, b, GAUSS)
     case(15)
-        call Solvers_Mod_Lu(grid, Ad, x, b, DOOLITTLE)
+        call Solvers_Mod_Lu(Grid, Ad, x, b, DOOLITTLE)
 
     case(21)
-        call Solvers_Mod_Incomplete_Cholesky(grid, As, x, b, f_in)
+        call Solvers_Mod_Incomplete_Cholesky(Grid, As, x, b, f_in)
     case(22)
-        call Solvers_Mod_Incomplete_Ldlt(grid, As, x, b, f_in)
+        call Solvers_Mod_Incomplete_Ldlt(Grid, As, x, b, f_in)
     case(23)
-        call Solvers_Mod_Incomplete_Lu(grid, As, x, b, f_in, GAUSS)
+        call Solvers_Mod_Incomplete_Lu(Grid, As, x, b, f_in, GAUSS)
     case(24)
-        call Solvers_Mod_Incomplete_Lu(grid, As, x, b, f_in, DOOLITTLE)
+        call Solvers_Mod_Incomplete_Lu(Grid, As, x, b, f_in, DOOLITTLE)
     case(25)
-        call Solvers_Mod_Incomplete_Ldlt_From_Tflows(grid, As, x, b)
+        call Solvers_Mod_Incomplete_Ldlt_From_Tflows(Grid, As, x, b)
 
     case(31)
-        call Solvers_Mod_Cg_No_Prec(grid, As, x, b, n_iter, res)
+        call Solvers_Mod_Cg_No_Prec(Grid, As, x, b, n_iter, res)
     case(32)
-        call Solvers_Mod_Cg_Diag_Prec(grid, As, x, b, n_iter, res)
+        call Solvers_Mod_Cg_Diag_Prec(Grid, As, x, b, n_iter, res)
     case(33)
-        call Solvers_Mod_Cg_Tflows_Prec(grid, As, x, b, n_iter, res)
+        call Solvers_Mod_Cg_Tflows_Prec(Grid, As, x, b, n_iter, res)
     case(34)
-        call Solvers_Mod_Cg_Ldlt_Prec(grid, As, x, b, n_iter, res, f_in)
+        call Solvers_Mod_Cg_Ldlt_Prec(Grid, As, x, b, n_iter, res, f_in)
     case(35)
-        call Solvers_Mod_Cg_Cholesky_Prec(grid, As, x, b, n_iter, res, f_in)
+        call Solvers_Mod_Cg_Cholesky_Prec(Grid, As, x, b, n_iter, res, f_in)
 
     case(41)
         call Demo_Mod_Compress_Decompress
     case(42)
-        call Demo_Mod_Fill_In(f_in, grid)
+        call Demo_Mod_Fill_In(f_in, Grid)
 
     case(51)
         print *, "# Enter the desired resolution: "
-        read *, grid % nx, grid % ny, grid % nz
+        read *, Grid % nx, Grid % ny, Grid % nz
     case(52)
         print *, "# Enter the desired fill-in level: "
         read *, f_in

Discretize_Mod/On_Sparse_Matrix.f90

@@ -1,10 +1,10 @@
 !==============================================================================!
-  subroutine On_Sparse_Matrix(Discrete, grid, A, x, b)
+  subroutine On_Sparse_Matrix(Discrete, Grid, A, x, b)
 !------------------------------------------------------------------------------!
   implicit none
 !---------------------------------[Arguments]----------------------------------!
   class(Discretize_Type) :: Discrete
-  type(Grid_Type)        :: grid
+  type(Grid_Type)        :: Grid
   type(Sparse_Type)      :: A
   real, allocatable      :: x(:)
   real, allocatable      :: b(:)
@@ -13,7 +13,7 @@ subroutine On_Sparse_Matrix(Discrete, grid, A, x, b)
 !==============================================================================!
 
   ! Create sparse system matrix
-  call A % Sparse_Create(grid)
+  call A % Sparse_Create(Grid, singular=.false.)
 
   call IO % Plot_Sparse ("sparse_a",  A)
   call IO % Print_Sparse("Sparse A:", A)

Grid_Mod.f90

@@ -7,7 +7,7 @@ module Grid_Mod
 !------------------------------------------------------------------------------!
   implicit none
 !------------------------------------------------------------------------------!
-!   Defines grid (which is ludicrously simple)                                 !
+!   Defines a (ludicrously simple) grid                                        !
 !==============================================================================!
 
   !---------------!
@@ -16,8 +16,29 @@ module Grid_Mod
   !               !
   !---------------!
   type Grid_Type
-    real    :: lx, ly, lz  ! domain sizes in x, y and z directions
-    integer :: nx, ny, nz  ! domain resolutions in x, y and z directions
+
+    integer :: n_cells
+    integer :: n_bnd_cells
+    integer :: n_faces
+
+    integer, allocatable :: faces_c(:,:)
+
+    integer :: nx, ny, nz  ! domain resolution in x, y and z direction
+    real    :: lx, ly, lz  ! domain size in x, y and z direction
+    real    :: dx, dy, dz  ! cell size in x, y and z direction
+
+    contains
+      procedure :: Create_Grid
+      procedure :: Cell_Number
+      procedure :: Cells_I_J_K
+      procedure :: Save_Vtk_Debug
+
   end type
 
+  contains
+#   include "Grid_Mod/Create_Grid.f90"
+#   include "Grid_Mod/Cell_Number.f90"
+#   include "Grid_Mod/Cells_I_J_K.f90"
+#   include "Grid_Mod/Save_Vtk_Debug.f90"
+
   end module

Grid_Mod/Cell_Number.f90

@@ -0,0 +1,21 @@
+!==============================================================================!
+  integer function Cell_Number(Grid, i, j, k)
+!------------------------------------------------------------------------------!
+!>  Returns the cell number based on i, j and k indices in Cartesian grid
+!------------------------------------------------------------------------------!
+  implicit none
+!---------------------------------[Arguments]----------------------------------!
+  class(Grid_Type)    :: Grid
+  integer, intent(in) :: i, j, k
+!-----------------------------------[Locals]-----------------------------------!
+  integer :: ni, nj, nk
+!==============================================================================!
+
+  ni = Grid % nx
+  nj = Grid % ny
+  nk = Grid % nz
+
+  Cell_Number = (k-1)*ni*nj + (j-1)*ni + i
+
+  end function
+

Grid_Mod/Cells_I_J_K.f90

@@ -0,0 +1,17 @@
+!==============================================================================!
+  subroutine Cells_I_J_K(Grid, c, i, j, k)
+!------------------------------------------------------------------------------!
+!>  Returns the cell's indices i, j and k based on the cell number.
+!------------------------------------------------------------------------------!
+  implicit none
+!---------------------------------[Arguments]----------------------------------!
+  class(Grid_Type)     :: Grid
+  integer, intent(in)  :: c
+  integer, intent(out) :: i, j, k
+!==============================================================================!
+
+  i = mod(c-1, Grid % nx) + 1
+  j = mod(c-i, Grid % nx * Grid % ny) / Grid % nx + 1
+  k = (c - i - (j-1) * Grid % nx) / (Grid % nx * Grid % ny) + 1
+
+  end subroutine