fix for rotation around 90 deg

more tests

python/damask/_geom.py

@@ -590,19 +590,23 @@ class Geom:
         dtype = float if np.isnan(fill) or int(fill) != fill or self.microstructure.dtype==np.float else int
 
         Eulers = R.as_Eulers(degrees=True)
-
         microstructure_in = self.get_microstructure()
+
         # These rotations are always applied in the reference coordinate system, i.e. (z,x,z) not (z,x',z'')
         # see https://www.cs.utexas.edu/~theshark/courses/cs354/lectures/cs354-14.pdf
-        for angle,axis in zip(Eulers[::-1], [(0,1),(1,2),(0,1)]):
-            microstructure_out = ndimage.rotate(microstructure_in,angle,axis,order=0,
+        for angle,axes in zip(Eulers[::-1], [(0,1),(1,2),(0,1)]):
+            microstructure_out = ndimage.rotate(microstructure_in,angle,axes,order=0,
                                                 prefilter=False,output=dtype,cval=fill)
-            microstructure_in = microstructure_out
+            if np.prod(microstructure_in.shape) == np.prod(microstructure_out.shape):
+                # avoid scipy interpolation errors for rotations close to multiples of 90°
+                microstructure_in = np.rot90(microstructure_in,k=np.rint(angle/90.).astype(int),axes=axes)
+            else:
+                microstructure_in = microstructure_out
 
-        origin = self.origin-(np.asarray(microstructure_out.shape)-self.grid)//2 * self.size/self.grid
+        origin = self.origin-(np.asarray(microstructure_in.shape)-self.grid)*.5 * self.size/self.grid
 
         #self.add_comments('geom.py:renumber v{}'.format(version)
-        return self.update(microstructure_out,origin=origin,rescale=True)
+        return self.update(microstructure_in,origin=origin,rescale=True)
 
 
     def canvas(self,grid=None,offset=None,fill=None):

python/tests/reference/Geom/rotate_Eulers=0.0-32.0-240.0.geom

@@ -0,0 +1,85 @@
+4 header
+grid   a 8 b 10 c 8
+size   x 8e-06 y 1e-05 z 8e-06
+origin x 0.0 y -2.5e-06 z -2e-06
+homogenization 1
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 30 14 18 42 42 42
+42 42 29 13 17 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 31 15 19 42 42 42
+42  6 10  2  2 42 42 42
+42  1  2  2  2  2 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 32 16 20 42 42 42
+42  7 11  2  2 42 42 42
+42  7 11  2  2 42 42 42
+42  2  2  2  2  2 42 42
+42 42  2  2 31 35 42 42
+42 42 22 26 10 14  1 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42  4  2  2  2  2 42 42
+42 42  2  2 32 36 42 42
+42 42 24 28 12 16  1 42
+42 42 42  7  7  1  1 42
+42 42 42  7  7  1  1 42
+42 42 42  1  1  1 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 25 29 13 17  1 42
+42 42 42  8  8  1  1 42
+42 42 42  1  1  1 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42  1  1  1 42 42
+42 42 42  1  1  1 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42

python/tests/reference/Geom/rotate_Eulers=32.0-68.0-21.0.geom

@@ -0,0 +1,104 @@
+4 header
+grid   a 11 b 11 c 9
+size   x 1.1e-05 y 1.1000000000000001e-05 z 9e-06
+origin x -1.5e-06 y -3.0000000000000005e-06 z -2.4999999999999998e-06
+homogenization 1
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42  1  2 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42  1 42 42 42 42 42 42
+42 42 42 42  1  5 42 42 42 42 42
+42 42 42  1  7  4 42 42 42 42 42
+42 42 42 42 42 27 42 42 42 42 42
+42 42 42 42 42 42  2 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42  1  1 42 42 42 42 42 42
+42 42 42  1  1  9 29 42 42 42 42
+42 42  1  1 11  8 28  2 42 42 42
+42 42 42  1 10 31  2 42 42 42 42
+42 42 42 42 30  2  2  2 42 42 42
+42 42 42 42 42 42  2  1 42 42 42
+42 42 42 42 42 42 42  1 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42  1 42 42 42 42 42 42 42
+42 42 42  1  1 42 42 42 42 42 42
+42 42  1 16 36 12 32 42 42 42 42
+42 42 42 15 35  2  2  2 42 42 42
+42 42 42 42  2  2  2 11  3 42 42
+42 42 42 42 42 42 10  6 42 42 42
+42 42 42 42 42 42 42  6 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42  1 42 42 42 42 42 42 42
+42 42 42  1 17 42 42 42 42 42 42
+42 42 42 16 36 37  2 42 42 42 42
+42 42 42 42 39  2  2 12 42 42 42
+42 42 42 38  2  2  2 11  8 42 42
+42 42 42 42  2  2 14 30 42 42 42
+42 42 42 42 42 42 13 30 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 40 42 42 42 42 42 42
+42 42 42 42 42  2 42 42 42 42 42
+42 42 42 42 42  2  2 15 42 42 42
+42 42 42 42 42  2 18 42 42 42 42
+42 42 42 42 42 42 17 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42  2 20 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42
+42 42 42 42 42 42 42 42 42 42 42

python/tests/test_Geom.py

@@ -5,12 +5,13 @@ import pytest
 import numpy as np
 
 from damask import Geom
+from damask import Rotation
 
 
 def geom_equal(a,b):
     return np.all(a.get_microstructure() == b.get_microstructure()) and \
-           np.all(a.get_size()           == b.get_size())           and \
-           np.all(a.get_grid()           == b.get_grid())
+           np.all(a.get_grid()           == b.get_grid()) and \
+           np.allclose(a.get_size(), b.get_size())
 
 @pytest.fixture
 def default():
@@ -98,6 +99,46 @@ class TestGeom:
         if update: modified.to_file(reference)
         assert geom_equal(modified,Geom.from_file(reference))
 
+    def test_renumber(self,default):
+        modified = copy.deepcopy(default)
+        microstructure = modified.get_microstructure()
+        for m in np.unique(microstructure):
+            microstructure[microstructure==m] = microstructure.max() + np.random.randint(1,30)
+        modified.update(microstructure)
+        assert not geom_equal(modified,default)
+        modified.renumber()
+        assert geom_equal(modified,default)
+
+    def test_substitute(self,default):
+        modified = copy.deepcopy(default)
+        microstructure = modified.get_microstructure()
+        offset = np.random.randint(1,500)
+        microstructure+=offset
+        modified.update(microstructure)
+        assert not geom_equal(modified,default)
+        modified.substitute(np.arange(default.microstructure.max())+1+offset,
+                            np.arange(default.microstructure.max())+1)
+        assert geom_equal(modified,default)
+
+    @pytest.mark.parametrize('axis_angle',[np.array([1,0,0,86.7]), np.array([0,1,0,90.4]), np.array([0,0,1,90]),
+                                           np.array([1,0,0,175]),np.array([0,-1,0,178]),np.array([0,0,1,180])])
+    def test_rotate90(self,default,axis_angle):
+        modified = copy.deepcopy(default)
+        for i in range(np.rint(360/axis_angle[3]).astype(int)):
+            modified.rotate(Rotation.from_axis_angle(axis_angle,degrees=True))
+        assert geom_equal(modified,default)
+
+    @pytest.mark.parametrize('Eulers',[[32.0,68.0,21.0],
+                                       [0.0,32.0,240.0]])
+    def test_rotate(self,default,update,reference_dir,Eulers):
+        modified = copy.deepcopy(default)
+        modified.rotate(Rotation.from_Eulers(Eulers,degrees=True))
+        tag = 'Eulers={}-{}-{}'.format(*Eulers)
+        reference = os.path.join(reference_dir,'rotate_{}.geom'.format(tag))
+        if update: modified.to_file(reference)
+        assert geom_equal(modified,Geom.from_file(reference))
+
+
     @pytest.mark.parametrize('periodic',[True,False])
     def test_tessellation_approaches(self,periodic):
         grid   = np.random.randint(10,20,3)