Search code examples
fortran

How to avoid code duplication between 2D and 3D versions?


I am writing two programs (a 2D and a 3D version). Nowadays, they are becoming huge (5000 code lines). There are many duplications and it is becoming tedious and difficult to maintain. So I decided to have only one source code that it could be compiled to the 2D or the 3D version.

To sum up, my problem, I wrote very simple examples (2D and 3D) but very representative of my problem.

2D version :
module two_dim_m

  implicit none
  
  type :: int_2dim_t
    integer :: x
    integer :: y
  end type int_2dim_t

end module two_dim_m

program main
  
  use two_dim_m
  
  implicit none
  type(int_2dim_t) :: ii_curr, full_dim
  
  full_dim%x = 5
  full_dim%y = 3
  
  associate (ii_x => ii_curr%x, ii_y => ii_curr%y)
  do ii_x = 1, full_dim%x
    do ii_y = 1, full_dim%y
      call my_print(ii_curr)
    enddo
  enddo
  end associate
  
contains

  subroutine my_print(ii)
    type(int_2dim_t), intent(in) :: ii
    
    write (*,'(i6,3x,i6,3x,f12.5)') ii%x,ii%y,sqrt(real(ii%x**2+ii%y**2))
    
  end subroutine my_print
    
end program main

3D version :

module three_dim_m

  implicit none
  
  type :: int_3dim_t
    integer :: x
    integer :: y
    integer :: z
  end type int_3dim_t

end module three_dim_m

program main
  
  use three_dim_m
  
  implicit none
  type(int_3dim_t) :: ii_curr, full_dim
  
  full_dim%x = 4
  full_dim%y = 3
  full_dim%z = 2
  
  associate (ii_x => ii_curr%x, ii_y => ii_curr%y, ii_z => ii_curr%z)
  do ii_x = 1, full_dim%x
    do ii_y = 1, full_dim%y
      do ii_z = 1, full_dim%z
        call my_print(ii_curr)
      enddo
    enddo
  enddo
  end associate

contains

  subroutine my_print(ii)
    type(int_3dim_t), intent(in) :: ii
    
    write (*,'(3(i6,3x),f12.5)') ii%x,ii%y,ii%z,sqrt(real(ii%x**2+ii%y**2+ii%z**2))
    
  end subroutine my_print
  
end program main

As, you could remark many lines are duplicated. Thus, I decided to use the C preprocessor to try to have only one code version without many line duplication (perhaps, not the best idea). After many tries (most of the ideas comes from stackoverflow answers), I wrote the following code :

! Just change dim_size to compile 2D or 3D version
#define cpp_dim_size 3

#if cpp_dim_size == 2
#define cpp_addDim(param) param%x, param%y
#define cpp_init_nLoop(currL,endL) \
  associate (ii_x => currL%x, ii_y => currL%y) cpp_nl \
  do ii_x = 1, endL%x cpp_nl \
    do ii_y = 1, endL%y
#define cpp_end_nLoop \
    enddo cpp_nl \
  enddo cpp_nl \
  end associate
  
#else if cpp_dim_size == 3
#define cpp_addDim(param) param%x, param%y, param%z
#define cpp_init_nLoop(currL,endL) \
  associate (ii_x => currL%x, ii_y => currL%y, ii_z => currL%z) cpp_nl \
  do ii_x = 1, endL%x cpp_nl \
    do ii_y = 1, endL%y cpp_nl \
      do ii_z = 1, endL%z
#define cpp_end_nLoop \
      enddo cpp_nl \
    enddo cpp_nl \
  enddo cpp_nl \
  end associate      
   
#endif

module n_dim_m

  implicit none
  
  type :: int_ndim_t
    integer :: x
    integer :: y
#if cpp_dim_size == 3
    integer :: z
#endif
  end type int_ndim_t

end module n_dim_m

program main
  
  use n_dim_m
  
  implicit none
  type(int_ndim_t) :: ii_curr, full_dim
  
  full_dim%x = 4
  full_dim%y = 3
#if cpp_dim_size == 3
    full_dim%z = 2
#endif

  cpp_init_nLoop(ii_curr,full_dim)
    call my_print(ii_curr)
  cpp_end_nLoop

contains

  subroutine my_print(ii)
    type(int_ndim_t), intent(in) :: ii
    character (len=32) :: fmt_string
    
    write (fmt_string,'(a,i1.1,a)') '(', cpp_dim_size, '(i6,3x),f12.5)'
    
    write (*, fmt= fmt_string) cpp_addDim(ii), &
#if cpp_dim_size == 2
            sqrt(real(ii%x**2+ii%y**2))
#else if cpp_dim_size == 3
            sqrt(real(ii%x**2+ii%y**2+ii%z**2))
#endif

  end subroutine my_print

end program main

To compile it :

gfortran  -cpp -P  -E fortran_cpp2.f90 | sed -e 's/cpp_nl/\n/g' >aaa.f90 ; gfortran aaa.f90

(Because of length line max size in fortran, sed is used to change cpp_nl to a line break)

It works, but I am not sure that it is the best solution. I tried to use OOP and polymorphism but I could not avoid line code duplication.


Solution

  • If you just need 2D and 3D code (i.e. you don't need to go higher than 3D), then the simplest solution is probably to store everything using 3D arrays, but to have the final dimension of these arrays be of size 1 for the 2D case.

    Something like

    module m
      implicit none
      
      type :: int_t
        integer, allocatable :: point(3)
      end type int_t
    end module m
    
    program main
      use m
      implicit none
      
      type(int_t) :: ii_curr, full_dim
      
      ! You might want a 2D constructor, to save you having to write the `1`s.
      full_dim%point = [5, 3, 1]
      
      associate (ii_x => ii_curr%point(1), &
                &ii_y => ii_curr%point(2), &
                &ii_z => ii_curr%point(3))
      do ii_x = 1, full_dim%x
        do ii_y = 1, full_dim%y
          do ii_z = 1, full_dim%z
            call my_print(ii_curr)
          enddo
        enddo
      enddo
      end associate
    contains
    
      subroutine my_print(ii)
        type(int_t), intent(in) :: ii
        
        write (*,'(i6,3x,i6,3x,i6,3x,f12.5)') ii%x,ii%y,ii%z,sqrt(real(ii%x**2+ii%y**2+ii%z**2))
      end subroutine my_print
        
    end program main
    

    The triple loop there would be in the wrong order if z was always going to be of size 1, so you'd probably want to transpose it.

    You probably also want to keep a flag around telling you what your dimension actually is, to avoid the performance costs of the third dimension (but only where those performance costs actually appear - no point complicating your code for optimisations which don't speed your code up).