test_grid2d.f

Go to the documentation of this file.

C> @file tests/test_grid2d.f
C> Unit tests of subroutines from grid2D/grid2d.f
C GRIDMAN, grid managment library. Author: Vladislav Kotov, v.kotov@fz-juelich.de

!    Copyright (c) 2017  Forschungszentrum Juelich GmbH
!                        Vladislav Kotov 
!
!    This file is part of GRIDMAN.
!
!    GRIDMAN is free software: you can redistribute it and/or modify
!    it under the terms of the GNU General Public License as published by
!    the Free Software Foundation, either version 3 of the License, or
!    (at your option) any later version.
!
!    GRIDMAN is distributed in the hope that it will be useful,
!    but WITHOUT ANY WARRANTY; without even the implied warranty of
!    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
!    GNU General Public License for more details.
!
!    You should have received a copy of the GNU General Public License
!    along with GRIDMAN.  If not, see <http://www.gnu.org/licenses/>.

      PROGRAM test_grid2d
      USE gridman
      USE gridman_lib
      IMPLICIT NONE

      INTRINSIC any,sum,dot_product,sqrt

      TYPE(gridman_grid) :: grid2d
      INTEGER :: ierr,res
      REAL(GRIDMAN_DP),ALLOCATABLE :: ledges(:),sedges(:),
     r                                vol(:),xcn(:,:),
     r                                scell(:),vn(:,:)
      REAL(GRIDMAN_DP) :: ltot,stot,vtot,xtmp,xc,yc
      LOGICAL,ALLOCATABLE :: isconvex(:)

      gridman_dbg=.true.
      gridman_dbg=.false.

      CALL grid_example1(grid2d,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"
      CALL gridman_grid2d_check(grid2d,res,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"
      IF(res.NE.0) THEN
       WRITE(*,*) "ERROR in GRIDMAN_GRID2D_CHECK"
       WRITE(*,*) "Expected value 0, RES ",res
       stop "TEST_GRID2D TERMINATED"
      END IF

C     TEST _GRID2D_CHECK
      CALL test_check(grid2d)

      gridman_check=.true.

C
C     NO SEPARATE TEST FOR GRIDMAN_GRID2D_CHAINS   
C     - TESTED INSIDE OTHER SUBROUTINES
C

C     EDGE LENGTHS 
      ALLOCATE(ledges(grid2d%NEDGES))
      CALL gridman_grid2d_lengths(grid2d,ledges,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"
      ltot=sum(ledges)
      IF(abs(ltot-58.155374460927121_gridman_dp).GT.gridman_tol) THEN
       WRITE(*,*) "ERORR in GRIDMAN_GRID2D_LENGTHS ",
     w            " wrong total edge length"
       WRITE(*,*) " LTOT ",ltot      
       stop "TEST_GRID2D TERMINATED"
      END IF
      DEALLOCATE(ledges)

C     CELL CROSS SECTIONS
      ALLOCATE(scell(grid2d%NCELLS))
      CALL gridman_grid2d_crossect(grid2d,scell,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"
      stot=sum(scell)
      IF(abs(stot-38.0).GT.gridman_tol) THEN
       WRITE(*,*) "ERORR in GRIDMAN_GRID2D_CELLAREAS",
     w            " wrong total poloidal area"
       WRITE(*,*) " STOT ",stot      
       stop "TEST_GRID2D TERMINATED"
      END IF

C     CELL CENTERS
      ALLOCATE(xcn(2,grid2d%NCELLS))
      CALL gridman_grid2d_center(grid2d,xcn,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"
cc      DO I=1,GRID2D%NCELLS
cc       WRITE(*,'(1X,I3,2F7.3)') I,XCN(1,I),XCN(2,I)
cc      END DO   
      xc=dot_product(xcn(1,:),scell)
      yc=dot_product(xcn(2,:),scell)
      stot=sum(scell)
      xc=xc/stot
      yc=yc/stot
       IF(abs(xc-7.0).GT.gridman_tol.OR.
     w    abs(yc+2.105263157894737e-01).GT.gridman_tol) THEN
       WRITE(*,*) "ERORR in GRIDMAN_GRID2D_CENTER",
     w            " wrong center of cells coordinates"
       WRITE(*,*) "XC, YC ",xc,yc
       stop "TEST_GRID2D TERMINATED"
      END IF
      DEALLOCATE(xcn)
      DEALLOCATE(scell)

C     AREAS 
      ALLOCATE(sedges(grid2d%NEDGES))
      CALL gridman_grid2d_cylareas(grid2d,sedges,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"
      stot=sum(sedges)
      IF(abs(stot-2557.8070316024337_gridman_dp).GT.1e-4) THEN
       WRITE(*,*) "ERORR in GRIDMAN_GRID2D_AREAS:",
     w            " wrong edge areas"
       WRITE(*,*) " STOT ",stot
       stop "TEST_GRID2D TERMINATED"
      END IF
      CALL gridman_grid2d_cylareas(grid2d,sedges,ierr,gridman_pi)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"
      stot=sum(sedges)
      IF(abs(stot-0.5*2557.8070316024337_gridman_dp).GT.1e-4) THEN
       WRITE(*,*) "ERORR in GRIDMAN_GRID2D_AREAS:",
     w            " wrong edge areas"
       WRITE(*,*) " STOT ",stot
       stop "TEST_GRID2D TERMINATED"
      END IF
      DEALLOCATE(sedges) 

C     VOLUMES 
      ALLOCATE(vol(grid2d%NCELLS))
      CALL gridman_grid2d_cylvolumes(grid2d,vol,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"
      vtot=sum(vol)
      IF(abs(vtot-1671.3273382186890_gridman_dp).GT.1e-4) THEN
       WRITE(*,*) "ERORR in GRIDMAN_GRID2D_VOLUMES",
     w            " wrong total volume"
       WRITE(*,*) " VTOT ",vtot
       stop "TEST_GRID2D TERMINATED"
      END IF 
      CALL gridman_grid2d_cylvolumes(grid2d,vol,ierr,0.5*gridman_pi)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"
      vtot=sum(vol)
      IF(abs(vtot-0.25*1671.3273382186890_gridman_dp).GT.1e-4) THEN
       WRITE(*,*) "ERORR in GRIDMAN_GRID2D_VOLUMES",
     w            " wrong total volume"
       WRITE(*,*) " VTOT ",vtot
       stop "TEST_GRID2D TERMINATED"
      END IF 
      DEALLOCATE(vol)
      
C     CELLS CONVEX OR NOT?
      ALLOCATE(isconvex(grid2d%NCELLS))
      CALL gridman_grid2d_isconvex(grid2d,isconvex,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"
cc      write(*,*) "isconvex ",isconvex
      IF(.NOT.all(isconvex)) THEN
       WRITE(*,*) "ERORR in GRIDMAN_GRID2D_ISCONVEX",
     w           " shows that not all cells are convex"       
       stop "TEST_GRID2D TERMINATED"
      END IF 

C     SHIFT ONE POINT TO CREATE A NON-CONVEX CELL
      xtmp=grid2d%X(2,10)
      grid2d%X(2,10)=-3.0
      CALL gridman_grid2d_isconvex(grid2d,isconvex,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"
cc      write(*,*) "isconvex ",isconvex
      IF(.NOT.all(isconvex(1:9))) THEN
       WRITE(*,*) "ERORR in GRIDMAN_GRID2D_ISCONVEX: ",
     w   " shows that not all cells from 1st to 9th are convex"       
       stop "TEST_GRID2D TERMINATED"
      END IF
      IF(isconvex(10)) THEN
       WRITE(*,*) "ERORR in GRIDMAN_GRID2D_ISCONVEX",
     w            " shows that 10th cell is convex"       
       stop "TEST_GRID2D TERMINATED"
      END IF
      grid2d%X(2,10)=xtmp
      DEALLOCATE(isconvex)

      write(*,*) "GRIDMAN_GRID2D_NORM "
C     EDGE NORMALS
      ALLOCATE(vn(2,grid2d%NEDGES))
      CALL gridman_grid2d_norm(grid2d,vn,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"      
cc      DO I=1,GRID2D%NEDGES
cc       WRITE(*,*) I,VN(1,I),VN(2,I)
cc      END DO  
      CALL test_normals(vn,grid2d%NEDGES)
      DEALLOCATE(vn)

      CALL gridman_grid_deallocate(grid2d,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"

      WRITE(*,*) "TEST_GRID2D COMPLETED"

      END PROGRAM test_grid2d

C
C
C
      SUBROUTINE test_normals(VN,N)
      USE gridman
      IMPLICIT NONE      
      INTRINSIC abs
      INTEGER(GRIDMAN_DP),INTENT(IN) :: n
       REAL(GRIDMAN_DP),INTENT(IN) :: vn(2,n)

      INTEGER :: i      

      i=1
      IF(abs(vn(1,i)).GT.0..OR.abs(vn(2,i)-1.0).GT.0.) GOTO 100
      i=2
      IF(abs(vn(1,i)).GT.0..OR.abs(vn(2,i)-1.0).GT.0.) GOTO 100
      i=3
      IF(abs(vn(1,i)+1.0).GT.0..OR.abs(vn(2,i)).GT.0.) GOTO 100
      i=5
      IF(abs(vn(1,i)-1.0).GT.0..OR.abs(vn(2,i)).GT.0.) GOTO 100
      i=8
      IF(abs(vn(1,i)+1.0).GT.0..OR.abs(vn(2,i)).GT.0.) GOTO 100
      i=10
      IF(abs(vn(1,i)-1.0).GT.0..OR.abs(vn(2,i)).GT.0.) GOTO 100
      i=15
      IF(abs(vn(1,i)).GT.0..OR.abs(vn(2,i)+1.0).GT.0.) GOTO 100
      i=16
      IF(abs(vn(1,i)).GT.0..OR.abs(vn(2,i)+1.0).GT.0.) GOTO 100
      i=21
      IF(abs(vn(1,i)+1.0).GT.0..OR.abs(vn(2,i)).GT.0.) GOTO 100
      i=22
      IF(abs(vn(1,i)-1.0).GT.0..OR.abs(vn(2,i)).GT.0.) GOTO 100
      i=23
      IF(abs(vn(1,i)).GT.0..OR.abs(vn(2,i)+1.0).GT.0.) GOTO 100
      i=13
      IF(abs(vn(1,i)+0.554700196225229).GT.gridman_tol.OR.
     f   abs(vn(2,i)-0.832050294337844).GT.gridman_tol) GOTO 100
      i=14
      IF(abs(vn(1,i)-0.554700196225229).GT.gridman_tol.OR.
     f   abs(vn(2,i)-0.832050294337844).GT.gridman_tol) GOTO 100

      RETURN

 100  WRITE(*,*) "ERORR in GRIDMAN_GRID2D_NORM ",
     w           "wrong direction of normal vector I ",i
      stop "TEST_GRID2D TERMINATED"

      END  SUBROUTINE test_normals

C
C    
C      
      SUBROUTINE test_check(GRID2D)
      USE gridman
      USE gridman_lib
      IMPLICIT NONE

      TYPE(gridman_grid) :: grid2d
      TYPE(gridman_grid) :: grid_copy
      INTEGER :: ierr,res
      INTEGER(GRIDMAN_SP) :: ip1,ip2
      
      CALL gridman_grid_copy(grid_copy,grid2d,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"

C     GRID TYPE
      grid2d%TYPE=0 
      CALL gridman_grid2d_check(grid2d,res,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"
      IF(res.NE.100) THEN
       WRITE(*,*) "ERROR in GRIDMAN_GRID2D_CHECK"
       WRITE(*,*) "Expected value 100, RES ",res
       stop "TEST_GRID2D TERMINATED"
      END IF
      CALL gridman_grid_copy(grid2d,grid_copy,ierr) !RESTORE OBJECT
      grid2d%PDIM=3
      CALL gridman_grid2d_check(grid2d,res,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"
      IF(res.NE.100) THEN
       WRITE(*,*) "ERROR in GRIDMAN_GRID2D_CHECK"
       WRITE(*,*) "Expected value 100, RES ",res
       stop "TEST_GRID2D TERMINATED"
      END IF
      CALL gridman_grid_copy(grid2d,grid_copy,ierr) !RESTORE OBJECT
      grid2d%EDIM=3
      CALL gridman_grid2d_check(grid2d,res,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"
      IF(res.NE.100) THEN
       WRITE(*,*) "ERROR in GRIDMAN_GRID2D_CHECK"
       WRITE(*,*) "Expected value 100, RES ",res
       stop "TEST_GRID2D TERMINATED"
      END IF
      CALL gridman_grid_copy(grid2d,grid_copy,ierr) !RESTORE OBJECT

C     _GRID_CHECK ERROR
      grid2d%EDIM=-1
      CALL gridman_grid2d_check(grid2d,res,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"
      IF(res.GT.100.OR.res.LT.1) THEN
       WRITE(*,*) "ERROR in GRIDMAN_GRID2D_CHECK"
       WRITE(*,*) "Expected value 0<R<100, RES ",res
       stop "TEST_GRID2D TERMINATED"
      END IF
      CALL gridman_grid_copy(grid2d,grid_copy,ierr) !RESTORE OBJECT

C     BROKEN CHAIN
      grid2d%POINTS(1,23)=12
      CALL gridman_grid2d_check(grid2d,res,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"
      IF(res.NE.102) THEN
       WRITE(*,*) "ERROR in GRIDMAN_GRID2D_CHECK"
       WRITE(*,*) "Expected value 102, RES ",res
       stop "TEST_GRID2D TERMINATED"
      END IF
      CALL gridman_grid_copy(grid2d,grid_copy,ierr) !RESTORE OBJECT

C     EDGE WITH ZERO LENGTH     
      ip1=grid2d%POINTS(1,23)
      ip2=grid2d%POINTS(2,23)
      grid2d%X(1,ip1)=grid2d%X(1,ip2)
      grid2d%X(2,ip1)=grid2d%X(2,ip2)
      CALL gridman_grid2d_check(grid2d,res,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"
      IF(res.NE.103) THEN
       WRITE(*,*) "ERROR in GRIDMAN_GRID2D_CHECK"
       WRITE(*,*) "Expected value 103, RES ",res
       stop "TEST_GRID2D TERMINATED"
      END IF
      CALL gridman_grid_copy(grid2d,grid_copy,ierr) !RESTORE OBJECT

C     INTERSECTING EDGES
      grid2d%POINTS(1,21)=11
      grid2d%POINTS(2,21)=16
      grid2d%POINTS(1,22)=15
      grid2d%POINTS(2,22)=12
      CALL gridman_grid2d_check(grid2d,res,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"
      IF(res.NE.111) THEN
       WRITE(*,*) "ERROR in GRIDMAN_GRID2D_CHECK"
       WRITE(*,*) "Expected value 111, RES ",res
       stop "TEST_GRID2D TERMINATED"
      END IF
      CALL gridman_grid_copy(grid2d,grid_copy,ierr) !RESTORE OBJECT

      CALL gridman_grid_deallocate(grid_copy,ierr)
      IF(ierr.NE.0) stop "TEST_GRID2D TERMINATED"

      END SUBROUTINE test_check