test_tria.f

Go to the documentation of this file.

C> @file tests/test_tria.f
C> Unit tests of subroutines from formats/tria.f and convert/tria.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_tria
       USE gridman
       USE gridman_lib
       IMPLICIT NONE    
       INTRINSIC sum

       CHARACTER(256) :: fnames(3)     
       INTEGER(GRIDMAN_SP) :: nrknot,ntrii,i,iy,icell
       INTEGER :: ierr 
       REAL(GRIDMAN_DP),ALLOCATABLE :: xytrian(:,:),xytrian2(:,:)
       INTEGER(GRIDMAN_SP),ALLOCATABLE :: necke(:,:),nchbar(:,:),
     i                                    nseite(:,:),itri(:,:),
     i                                    necke2(:,:),nchbar2(:,:),
     i                                    nseite2(:,:),itri2(:,:)
       TYPE(gridman_grid) :: tria_grid,tria_grid_copy,tria1,tria2,tria
       REAL(GRIDMAN_DP),ALLOCATABLE :: vol(:),vol1(:),vol2(:) 
       REAL(GRIDMAN_DP) :: v1,v2,v0,v
       LOGICAL,ALLOCATABLE :: lcadd(:)

       gridman_check=.true.
cc       GRIDMAN_DBG=.TRUE.

C      READ TRIANGULAR GRID
       fnames(1)='./input/fort.33';
       fnames(2)='./input/fort.34';
       fnames(3)='./input/fort.35';
       CALL gridman_tria_read(fnames,nrknot,ntrii,xytrian,
     s                        necke,nchbar,nseite,itri,ierr)
       IF(ierr.NE.0) stop "TEST_TRIA TERMINATED"
C      CONVERT TRIA INTO GRID OBJECT
       CALL gridman_tria2grid(tria_grid,nrknot,ntrii,xytrian,
     s                        necke,nchbar,nseite,itri,ierr)
       IF(ierr.NE.0) stop "TEST_TRIA TERMINATED" 
C      WRITE TRIANGULAR GRID
       fnames(1)='fort.33b';
       fnames(2)='fort.34b';
       fnames(3)='fort.35b';
       CALL gridman_tria_write(tria_grid,fnames,ierr)
       IF(ierr.NE.0) stop "TEST_TRIA TERMINATED"
       fnames(1)='';
       fnames(2)='';
       fnames(3)='';
       CALL gridman_tria_write(tria_grid,fnames,ierr)
       IF(ierr.NE.0) stop "TEST_TRIA TERMINATED"    
       CALL gridman_grid_write(tria_grid,'tria.grd',ierr)
       IF(ierr.NE.0) stop "TEST_TRIA TERMINATED"  

C      READ TRIANGULAR GRID WHICH WAS WRITTEN
       CALL gridman_tria_read(tria_grid_copy,fnames,ierr)
       IF(ierr.NE.0) stop "TEST_TRIA TERMINATED"

C      CHECK TRANSFORMATIONS OF ITRI
       ALLOCATE(xytrian2(2,tria_grid_copy%NPOINTS),
     a          necke2(3,tria_grid_copy%NCELLS),
     a          nchbar2(3,tria_grid_copy%NCELLS),
     a          nseite2(3,tria_grid_copy%NCELLS),
     a          itri2(5,tria_grid_copy%NCELLS)) 
       itri2=0
       CALL gridman_grid2tria(tria_grid_copy,xytrian2,necke2,nchbar2,
     c                        nseite2,itri2,ierr)
       IF(ierr.NE.0) stop "TEST_TRIA TERMINATED"
       DO icell=1,tria_grid_copy%NCELLS
        IF(itri(1,icell).NE.itri2(1,icell).OR.
     f     itri(2,icell).NE.itri2(2,icell)) THEN 
         WRITE(*,*) "ERROR: cell index mismatch"
         WRITE(*,*) " ICELL ",icell
         WRITE(*,*) " IX, IY ",itri(1,icell),itri(2,icell)
         WRITE(*,*) " IX2, IY2 ",itri2(1,icell),itri2(2,icell)
         stop "TEST_TRIA TERMINATED" 
        END IF
        IF(sum(itri(3:5,icell)).NE.sum(itri2(3:5,icell))) THEN
         WRITE(*,*) "ERROR: edge index mismatch"
         WRITE(*,*) " ICELL ",icell
         WRITE(*,*) " ITRI ",itri(3:5,icell)
         WRITE(*,*) " ITRI2 ",itri2(3:5,icell)
         stop "TEST_TRIA TERMINATED" 
        END IF
       END DO
       DEALLOCATE(xytrian2,necke2,nchbar2,nseite2,itri2)   

C      CHECK IF TRIA_GRID AND TRIA_GRID_COPY HAVE SAME VOLUMES
       ALLOCATE(vol1(tria_grid%NCELLS),vol2(tria_grid_copy%NCELLS))
       gridman_tol=1e-5
       CALL gridman_grid2d_cylvolumes(tria_grid,vol1,ierr)
       IF(ierr.NE.0) stop "TEST_TRIA TERMINATED" 
       CALL gridman_grid2d_cylvolumes(tria_grid_copy,vol2,ierr)
       IF(ierr.NE.0) stop "TEST_TRIA TERMINATED" 
       v1=sum(vol1)
       v2=sum(vol2)
       IF(abs(v1-v2).GT.gridman_tol*(abs(v1)+abs(v2))) THEN
        WRITE(*,*) "ERROR: grids have different volumes"
        WRITE(*,*) "V1, V2 ",v1,v2
        stop "TEST_TRIA TERMINATED" 
       END IF

       gridman_tol=1e-7
C      CONVERT INTO GRID OBJECT WHICH INCLUDES ONLY IY<11 CELLS
       ALLOCATE(lcadd(tria_grid%NCELLS))
       lcadd=.false.
       DO i=1,ntrii
        iy=itri(2,i)
        IF(iy.LT.11) lcadd(i)=.true.
       END DO        
       CALL gridman_grid_eliminate_cells(tria1,tria_grid,lcadd,ierr)
       IF(ierr.NE.0) stop "TEST_TRIA TERMINATED"

C      CONVERT INTO GRID OBJECT WHICH INCLUDES ONLY IY>=11 CELLS
       lcadd=.false.
       DO i=1,ntrii
        iy=itri(2,i)
        IF(iy.GE.11) lcadd(i)=.true.
       END DO        
       CALL gridman_grid_eliminate_cells(tria2,tria_grid,lcadd,ierr)
       IF(ierr.NE.0) stop "TEST_TRIA TERMINATED"

C      CHECK THAT SUM OF VOLUMES OF THOSE GRIDS
C      IS EQUAL TO THE VOLUME OF COMBINED GRID     
       DEALLOCATE(vol1,vol2)
       ALLOCATE(vol(tria_grid%NCELLS),
     a          vol1(tria1%NCELLS),vol2(tria2%NCELLS))
       gridman_tol=1e-5
       CALL gridman_grid2d_cylvolumes(tria_grid,vol,ierr)
       IF(ierr.NE.0) stop "TEST_TRIA TERMINATED" 
       CALL gridman_grid2d_cylvolumes(tria1,vol1,ierr)
       IF(ierr.NE.0) stop "TEST_TRIA TERMINATED" 
       CALL gridman_grid2d_cylvolumes(tria2,vol2,ierr)
       IF(ierr.NE.0) stop "TEST_TRIA TERMINATED" 
       v0=sum(vol)
       v1=sum(vol1)
       v2=sum(vol2)
       IF(abs(v0-v1-v2).GT.gridman_tol) THEN
        WRITE(*,*) "ERROR: sum of volumes is not equal",
     w             " to the volume of combined grid"
        WRITE(*,*) "V, V1+V2 ",v0,v1+v2
        stop "TEST_TRIA TERMINATED" 
       END IF

C     MERGE TWO PARTIAL GRIDS
      CALL gridman_grid2d_merge(tria,tria1,tria2,1e-5_gridman_dp,ierr)
      IF(ierr.NE.0) stop "TEST_MERGE TERMINATED"
      IF(tria%NCELLS.NE.tria_grid%NCELLS) THEN
       WRITE(*,*) "ERROR: wrong merging"
       stop "TEST_TRIA TERMINATED"
      END IF
      CALL gridman_grid2d_cylvolumes(tria,vol,ierr)
      IF(ierr.NE.0) stop "TEST_TRIA TERMINATED" 
      v=sum(vol)
      IF(abs(v0-v).GT.gridman_tol) THEN
       WRITE(*,*) "ERROR: mismatch between volumes ", 
     w            " of original and merged grids"
       WRITE(*,*) "V0, V ",v0,v
       stop "TEST_TRIA TERMINATED" 
      END IF

C     DEALLOCATE
      CALL gridman_grid_deallocate(tria_grid,ierr) 
      IF(ierr.NE.0) stop "TEST_TRIA TERMINATED"
      CALL gridman_grid_deallocate(tria_grid_copy,ierr) 
      IF(ierr.NE.0) stop "TEST_TRIA TERMINATED"
      CALL gridman_grid_deallocate(tria1,ierr) 
      IF(ierr.NE.0) stop "TEST_TRIA TERMINATED"
      CALL gridman_grid_deallocate(tria2,ierr) 
      IF(ierr.NE.0) stop "TEST_TRIA TERMINATED"
      CALL gridman_grid_deallocate(tria,ierr) 
      IF(ierr.NE.0) stop "TEST_TRIA TERMINATED"

       WRITE(*,*) "TEST_TRIA COMPLETED"
   
      END PROGRAM test_tria