polishing

.gitlab-ci.yml

@@ -3,6 +3,7 @@ stages:
   - prepare
   - python
   - compile
+  - setup
   - fortran
   - performance
   - deploy
@@ -161,8 +162,10 @@ compile_Marc:
     - master
     - release
 
+
+###################################################################################################
 setup_grid:
-  stage: compile
+  stage: setup
   script:
     - module load $IntelCompiler $MPI_Intel $PETSc_Intel
     - BUILD_DIR=$(mktemp -d)
@@ -174,7 +177,7 @@ setup_grid:
     - release
 
 setup_mesh:
-  stage: compile
+  stage: setup
   script:
     - module load $IntelCompiler $MPI_Intel $PETSc_Intel
     - BUILD_DIR=$(mktemp -d)

python/damask/_grid.py

@@ -1,8 +1,8 @@
+import os
 import copy
+import warnings
 import multiprocessing as mp
 from functools import partial
-import os
-import warnings
 
 import numpy as np
 import pandas as pd
@@ -20,7 +20,7 @@ class Grid:
     Geometry definition for grid solvers.
 
     Create and manipulate geometry definitions for storage as VTK
-    image data files ('.vit' extension). A grid contains the
+    image data files ('.vti' extension). A grid contains the
     material ID (referring to the entry in 'material.yaml') and
     the physical size.
     """
@@ -166,9 +166,7 @@ class Grid:
             Grid-based geometry from file.
 
         """
-        if str(fname).endswith('.vtr'):
+        v = VTK.load(fname if str(fname).endswith('.vti') else str(fname)+'.vti')
-            warnings.warn('Support for vtr files will be removed in DAMASK 3.1.0', DeprecationWarning,2)
-        v = VTK.load(fname)
         comments = v.get_comments()
         cells = np.array(v.vtk_data.GetDimensions())-1
         bbox  = np.array(v.vtk_data.GetBounds()).reshape(3,2).T

python/damask/_vtk.py

@@ -1,4 +1,5 @@
 import os
+import warnings
 import multiprocessing as mp
 from pathlib import Path
 
@@ -69,8 +70,6 @@ class VTK:
         """
         Create VTK of type vtk.vtkRectilinearGrid.
 
-        This is the common type for grid solver results.
-
         Parameters
         ----------
         grid : iterable of int, len (3)
@@ -86,6 +85,7 @@ class VTK:
             VTK-based geometry without nodal or cell data.
 
         """
+        warnings.warn('Support for vtr files will be removed in DAMASK 3.1.0', DeprecationWarning,2)
         vtk_data = vtk.vtkRectilinearGrid()
         vtk_data.SetDimensions(*(np.array(grid)+1))
         coord = [np_to_vtk(np.linspace(origin[i],origin[i]+size[i],grid[i]+1),deep=True) for i in [0,1,2]]
@@ -179,8 +179,8 @@ class VTK:
         Parameters
         ----------
         fname : str or pathlib.Path
-            Filename for reading. Valid extensions are .vtr, .vtu, .vtp, and .vtk.
+            Filename for reading. Valid extensions are .vti, .vtr, .vtu, .vtp, and .vtk.
-        dataset_type : {'vtkRectilinearGrid', 'vtkUnstructuredGrid', 'vtkPolyData'}, optional
+        dataset_type : {'vtkImageData', ''vtkRectilinearGrid', 'vtkUnstructuredGrid', 'vtkPolyData'}, optional
             Name of the vtk.vtkDataSet subclass when opening a .vtk file.
 
         Returns

src/grid/discretization_grid.f90

@@ -164,9 +164,9 @@ subroutine readVTI(grid,geomSize,origin,material)
   integer,     dimension(:), intent(out), allocatable :: &
     material
 
-  character(len=:), allocatable :: fileContent, dataType, headerType, temp
+  character(len=:), allocatable :: fileContent, dataType, headerType
-  logical :: inFile,inGrid,gotCoordinates,gotCellData,compressed
+  logical :: inFile,inImage,gotCellData,compressed
-  integer :: fileUnit, myStat, coord
+  integer :: fileUnit, myStat
   integer(pI64) :: &
     fileLength, &                                                                                   !< length of the geom file (in characters)
     startPos, endPos, &
@@ -186,39 +186,37 @@ subroutine readVTI(grid,geomSize,origin,material)
   close(fileUnit)
 
   inFile         = .false.
-  inGrid         = .false.
+  inImage        = .false.
-  gotCoordinates = .false.
   gotCelldata    = .false.
 
 !--------------------------------------------------------------------------------------------------
-! interpret XML file
+! parse XML file
   startPos = 1_pI64
   do while (startPos < len(fileContent,kind=pI64))
     endPos = startPos + index(fileContent(startPos:),IO_EOL,kind=pI64) - 2_pI64
     if (endPos < startPos) endPos = len(fileContent,kind=pI64)                                      ! end of file without new line
 
-    if(.not. inFile) then
+    if (.not. inFile) then
       if(index(fileContent(startPos:endPos),'<VTKFile',kind=pI64) /= 0_pI64) then
         inFile = .true.
-        if(.not. fileFormatOk(fileContent(startPos:endPos))) call IO_error(error_ID = 844, ext_msg='file format')
+        if (.not. fileFormatOk(fileContent(startPos:endPos))) call IO_error(error_ID = 844, ext_msg='file format')
         headerType = merge('UInt64','UInt32',getXMLValue(fileContent(startPos:endPos),'header_type')=='UInt64')
         compressed  = getXMLValue(fileContent(startPos:endPos),'compressor') == 'vtkZLibDataCompressor'
       endif
     else
-      if (.not. inGrid) then
+      if (.not. inImage) then
         if (index(fileContent(startPos:endPos),'<ImageData',kind=pI64) /= 0_pI64) then
-          inGrid = .true.
+          inImage = .true.
-          print*, 'x'
+          call cellsSizeOrigin(grid,geomSize,origin,fileContent(startPos:endPos))
-          call cellsSizeOrigin(fileContent(startPos:endPos))
         endif
       else
-        if(index(fileContent(startPos:endPos),'<CellData>',kind=pI64) /= 0_pI64) then
+        if (index(fileContent(startPos:endPos),'<CellData>',kind=pI64) /= 0_pI64) then
           gotCellData = .true.
           do while (index(fileContent(startPos:endPos),'</CellData>',kind=pI64) == 0_pI64)
-            if(index(fileContent(startPos:endPos),'<DataArray',kind=pI64) /= 0_pI64 .and. &
+            if (index(fileContent(startPos:endPos),'<DataArray',kind=pI64) /= 0_pI64 .and. &
                  getXMLValue(fileContent(startPos:endPos),'Name') == 'material' ) then
 
-              if(getXMLValue(fileContent(startPos:endPos),'format') /= 'binary') &
+              if (getXMLValue(fileContent(startPos:endPos),'format') /= 'binary') &
                 call IO_error(error_ID = 844, ext_msg='format (material)')
               dataType = getXMLValue(fileContent(startPos:endPos),'type')
 
@@ -235,7 +233,7 @@ subroutine readVTI(grid,geomSize,origin,material)
       endif
     endif
 
-    if(gotCellData .and. gotCoordinates) exit
+    if (gotCellData) exit
     startPos = endPos + 2_pI64
 
   end do
@@ -245,19 +243,21 @@ subroutine readVTI(grid,geomSize,origin,material)
   if(any(geomSize<=0))                call IO_error(error_ID = 844, ext_msg='size')
   if(any(grid<1))                     call IO_error(error_ID = 844, ext_msg='grid')
   material = material + 1
-  if(any(material<0))                 call IO_error(error_ID = 844, ext_msg='material ID < 0')
+  if(any(material<1))                 call IO_error(error_ID = 844, ext_msg='material ID < 0')
 
   contains
 
   !------------------------------------------------------------------------------------------------
   !> @brief determine size and origin from coordinates
   !------------------------------------------------------------------------------------------------
-  subroutine cellsSizeOrigin(header)
+  subroutine cellsSizeOrigin(c,s,o,header)
 
-    character(len=*), intent(in) :: header
+    integer, dimension(3),     intent(out) :: c
+    real(pReal), dimension(3), intent(out) :: s,o
+    character(len=*),          intent(in) :: header
 
     character(len=:), allocatable :: temp
-    real(pReal), dimension(:), allocatable :: coords,delta,origin
+    real(pReal), dimension(:), allocatable :: delta
     integer, dimension(:), allocatable :: stringPos
     integer :: i
 
@@ -265,15 +265,17 @@ subroutine readVTI(grid,geomSize,origin,material)
     if (getXMLValue(header,'Direction') /= '1 0 0 0 1 0 0 0 1') &
       call IO_error(error_ID = 844, ext_msg = 'coordinate order')
 
-    temp = getXMLValue(header,'Origin')
+    temp = getXMLValue(header,'WholeExtent')
-    origin = [(IO_floatValue(temp,IO_stringPos(temp),i),i=1,3)]
+    if (any([(IO_floatValue(temp,IO_stringPos(temp),i),i=1,5,2)] /= 0)) &
+      call IO_error(error_ID = 844, ext_msg = 'coordinate start')
+    c = [(IO_floatValue(temp,IO_stringPos(temp),i),i=2,6,2)]
+
     temp = getXMLValue(header,'Spacing')
     delta = [(IO_floatValue(temp,IO_stringPos(temp),i),i=1,3)]
-    temp = getXMLValue(header,'WholeExtent')
+    s = delta * real(c,pReal)
-    grid = [(IO_floatValue(temp,IO_stringPos(temp),i),i=2,6,2)]
 
-
+    temp = getXMLValue(header,'Origin')
-    geomSize = delta * real(grid,pReal)
+    o = [(IO_floatValue(temp,IO_stringPos(temp),i),i=1,3)]
 
   end subroutine