fixed tutorials, use of code formatter

Author: guglielmopadula

docs/source/_tutorials/tutorial-7-cffd.html

@@ -59,10 +59,14 @@ class BFFD(CFFD):
         >>> assert np.isclose(np.linalg.norm(bffd.fun(bffd.ffd(original_mesh_points[:-4]))-b),np.array([0.]))
         >>> new_mesh_points = bffd.ffd(original_mesh_points)
     '''
-    def __init__(self, n_control_points=None, fixval=None, weight_matrix=None, mask=None):
-        super().__init__(n_control_points,None,fixval,weight_matrix,mask)
+    def __init__(self,
+                 n_control_points=None,
+                 fixval=None,
+                 weight_matrix=None,
+                 mask=None):
+        super().__init__(n_control_points, None, fixval, weight_matrix, mask)
 
         def linfun(x):
             return np.mean(x.reshape(-1, 3), axis=0)
 
-        self.fun = linfun
+        self.fun = linfun

pygem/bffd.py

@@ -59,29 +59,34 @@ class CFFD(FFD):
         >>> assert np.isclose(np.linalg.norm(fun(cffd.ffd(original_mesh_points[:-4]))-b),np.array([0.]),atol=1e-06)
         >>> new_mesh_points = cffd.ffd(original_mesh_points)
     """
-    def __init__(self, n_control_points=None, fun=None, fixval=None, weight_matrix=None, mask=None ):
+    def __init__(self,
+                 n_control_points=None,
+                 fun=None,
+                 fixval=None,
+                 weight_matrix=None,
+                 mask=None):
         super().__init__(n_control_points)
 
         if mask is None:
-            self.mask=np.full((*self.n_control_points,3), True, dtype=bool)
+            self.mask = np.full((*self.n_control_points, 3), True, dtype=bool)
         else:
-            self.mask=mask
+            self.mask = mask
 
         if fixval is None:
-            self.fixval=np.array([1.])    
+            self.fixval = np.array([1.])
         else:
-            self.fixval=fixval
-        
+            self.fixval = fixval
+
         if fun is None:
-            self.fun=lambda x: self.fixval
-        
+            self.fun = lambda x: self.fixval
+
         else:
-            self.fun=fun
+            self.fun = fun
 
         if weight_matrix is None:
-            self.weight_matrix=np.eye(np.sum(self.mask.astype(int)))
+            self.weight_matrix = np.eye(np.sum(self.mask.astype(int)))
 
-    def adjust_control_points(self,src_pts):
+    def adjust_control_points(self, src_pts):
         '''
         Adjust the FFD control points such that F(ffd(src_pts))=c
             
@@ -91,11 +96,17 @@ def adjust_control_points(self,src_pts):
         '''
 
         saved_parameters = self._save_parameters()
-        indices=np.arange(np.prod(self.n_control_points)*3)[self.mask.reshape(-1)]
-        A, b = self._compute_linear_map(src_pts, saved_parameters.copy(),indices)
+        indices = np.arange(np.prod(self.n_control_points) *
+                            3)[self.mask.reshape(-1)]
+        A, b = self._compute_linear_map(src_pts, saved_parameters.copy(),
+                                        indices)
         d = A @ saved_parameters[indices] + b
-        invM=np.linalg.inv(self.weight_matrix)
-        deltax = np.linalg.multi_dot([invM , A.T , np.linalg.inv(np.linalg.multi_dot([A, invM, A.T])) , (self.fixval - d)])
+        invM = np.linalg.inv(self.weight_matrix)
+        deltax = np.linalg.multi_dot([
+            invM, A.T,
+            np.linalg.inv(np.linalg.multi_dot([A, invM, A.T])),
+            (self.fixval - d)
+        ])
         saved_parameters[indices] = saved_parameters[indices] + deltax
         self._load_parameters(saved_parameters)
 
@@ -134,14 +145,15 @@ def _load_parameters(self, tmp):
         self.array_mu_y = tmp[:, :, :, 1]
         self.array_mu_z = tmp[:, :, :, 2]
 
-    def read_parameters(self,filename='parameters.prm'):
+    def read_parameters(self, filename='parameters.prm'):
         super().read_parameters(filename)
-        self.mask=np.full((*self.n_control_points,3), True, dtype=bool)
-        self.weight_matrix=np.eye(np.sum(self.mask.astype(int)))
+        self.mask = np.full((*self.n_control_points, 3), True, dtype=bool)
+        self.weight_matrix = np.eye(np.sum(self.mask.astype(int)))
+
 
 # I see that a similar function already exists in pygem.utils, but it does not work for inputs and outputs of different dimensions
 
-    def _compute_linear_map(self, src_pts, saved_parameters,indices):
+    def _compute_linear_map(self, src_pts, saved_parameters, indices):
         '''
         Computes the coefficient and the intercept of the linear map from the control points to the output.
         
@@ -171,4 +183,4 @@ def _compute_linear_map(self, src_pts, saved_parameters,indices):
                               rcond=None)  #computation of the linear map
         A = sol[0].T[:, :-1]  #coefficient
         b = sol[0].T[:, -1]  #intercept
-        return A, b
+        return A, b

pygem/cffd.py

@@ -65,8 +65,13 @@ class VFFD(CFFD):
         >>> assert np.isclose(np.linalg.norm(vffd.fun(new_mesh_points)-b),np.array([0.]),atol=1e-07)
 
     '''
-    def __init__(self, triangles, n_control_points=None, fixval=None, weight_matrix=None, mask=None ):
-        super().__init__(n_control_points,None,fixval,weight_matrix,mask)
+    def __init__(self,
+                 triangles,
+                 n_control_points=None,
+                 fixval=None,
+                 weight_matrix=None,
+                 mask=None):
+        super().__init__(n_control_points, None, fixval, weight_matrix, mask)
         self.triangles = triangles
         self.vweight = [1 / 3, 1 / 3, 1 / 3]
 
@@ -79,12 +84,12 @@ def volume(x):
 
     def adjust_control_points(self, src_pts):
         self.vweight = np.abs(self.vweight) / np.sum(np.abs(self.vweight))
-        mask_bak=self.mask.copy()        
+        mask_bak = self.mask.copy()
         diffvolume = self.fixval - self.fun(self.ffd(src_pts))
         for i in range(3):
-            self.mask=np.full((*self.n_control_points,3), False, dtype=bool)
-            self.mask[:,:,:,i]=mask_bak[:,:,:,i].copy()
+            self.mask = np.full((*self.n_control_points, 3), False, dtype=bool)
+            self.mask[:, :, :, i] = mask_bak[:, :, :, i].copy()
             self.weight_matrix = np.eye(np.sum(self.mask.astype(int)))
             self.fixval = self.fun(
                 self.ffd(src_pts)) + self.vweight[i] * (diffvolume)
-            super().adjust_control_points(src_pts)
+            super().adjust_control_points(src_pts)

pygem/vffd.py

@@ -16,10 +16,10 @@ def test_nothing_happens_0(self):
         def fun(x):
             x = x.reshape(-1)
             return A @ x
+
         new_mesh_points = cffd(original_mesh_points)
         assert np.array_equal(original_mesh_points, new_mesh_points)
 
-
     def test_nothing_happens(self):
         np.random.seed(0)
         cffd = CFFD()
@@ -57,7 +57,8 @@ def fun(x):
         cffd.adjust_control_points(original_mesh_points)
         new_mesh_points = cffd.ffd(original_mesh_points)
         assert np.isclose(np.linalg.norm(fun(new_mesh_points) - b),
-                          np.array([0.0]),atol=1e-7)
+                          np.array([0.0]),
+                          atol=1e-7)
 
     def test_interpolation(self):
         cffd = CFFD()
@@ -73,8 +74,10 @@ def fun(x):
         cffd.fixval = b
         cffd.fun = fun
         save_par = cffd._save_parameters()
-        indices=np.arange(np.prod(cffd.n_control_points)*3)[cffd.mask.reshape(-1)]
-        C, d = cffd._compute_linear_map(original_mesh_points, save_par.copy(),indices)
+        indices = np.arange(np.prod(cffd.n_control_points) *
+                            3)[cffd.mask.reshape(-1)]
+        C, d = cffd._compute_linear_map(original_mesh_points, save_par.copy(),
+                                        indices)
         for i in range(2 * len(indices)):
             tmp = np.random.rand(len(indices))
             save_par[indices] = tmp

tests/test_cffd.py

@@ -19,32 +19,29 @@ def custom_linear_constraint(x):
 print("The custom linear function on the non deformed points is", custom_linear_constraint(x))
 print("The custom linear function on the classic FFD deformed points is", custom_linear_constraint(x_def))
 from pygem.cffd import CFFD
-ffd=CFFD([3,3,1])
+ffd=CFFD([3,3,1],custom_linear_constraint,np.array([1.]))
 np.random.seed(0)
 ffd.array_mu_x=ffd.array_mu_x+0.5*np.random.rand(*ffd.array_mu_x.shape)
 ffd.array_mu_y=ffd.array_mu_x+0.5*np.random.rand(*ffd.array_mu_x.shape)
-ffd.linconstraint=custom_linear_constraint
-ffd.valconstraint=np.array([1.])
-ffd.indices=np.arange(3*3*3*1).reshape(3,3,1,3)
-ffd.indices=ffd.indices[:,:,:,:-1] #removing z indices
-ffd.indices=ffd.indices.reshape(-1).tolist()
-ffd.M=np.eye(len(ffd.indices))#no weighting
+ffd.mask[:,:,:,-1]=False #removing z masks
+ffd.weight_matrix=np.eye(np.sum(ffd.mask))#no weighting
+ffd.adjust_control_points(x)
 x_def=ffd(x)
 plt.plot(x_def[:,0],x_def[:,1],'o')
 print("The custom linear function on the constrained FFD deformed points is", custom_linear_constraint(x_def))
 from pygem.bffd import BFFD
 def mesh_points(num_pts = 2000):
     indices = np.arange(0, num_pts, dtype=float) + 0.5
+
     phi = np.arccos(1 - 2*indices/num_pts)
     theta = np.pi * (1 + 5**0.5) * indices
+
     return np.array([np.cos(theta) * np.sin(phi), np.sin(theta) * np.sin(phi), np.cos(phi)]).T
 mesh = mesh_points()
-ffd = BFFD([2, 2, 2])
-ffd.valconstraint=np.array([0.,0.,0.])
+ffd = BFFD([2, 2, 2],np.array([0.,0.,0.]))
 ffd.array_mu_x[1, 1, 1] = 2
 ffd.array_mu_z[1, 1, 1] = 0
-ffd.indices=np.arange(2*2*2*3).tolist()
-ffd.M=np.eye(len(ffd.indices))
+ffd.adjust_control_points(mesh)
 mesh_def=ffd(mesh)
 print(np.mean(mesh_def,axis=0))
 ax = plt.figure(figsize=(8,8)).add_subplot(111, projection='3d')
@@ -65,19 +62,19 @@ def mesh_points(num_pts = 2000):
 ax.plot_trisurf(points[:,0], points[:,1], points[:,2], triangles=faces, cmap=plt.cm.Spectral)
 from pygem.vffd import VFFD
 vffd=VFFD(faces,[2,2,2])
-initvolume=vffd.linconstraint(points)
-vffd.valconstraint=np.array([initvolume])
-vffd.indices=np.arange(2*2*2*3).tolist()
-vffd.M=np.eye(len(vffd.indices))
+initvolume=vffd.fun(points)
+vffd.fixval=np.array([initvolume])
 vffd.vweight=np.array([0,1,0])
 np.random.seed(0)
 vffd.array_mu_x=vffd.array_mu_x+0.5*np.random.rand(2,2,2)
 vffd.array_mu_y=vffd.array_mu_y+0.5*np.random.rand(2,2,2)
 vffd.array_mu_z=vffd.array_mu_z+0.5*np.random.rand(2,2,2)
+vffd.adjust_control_points(points)
 mesh_def=vffd(points)
 mesh_def=mesh_def.reshape(points.shape)
 print("Percentage difference from the original mesh is ", np.linalg.norm(mesh_def-points)/np.linalg.norm(points)*100)
 fig = plt.figure(figsize=plt.figaspect(0.5))
 ax = fig.add_subplot(1, 2, 1, projection='3d')
 ax.plot_trisurf(mesh_def[:,0], mesh_def[:,1], mesh_def[:,2], triangles=faces, cmap=plt.cm.Spectral)
 
+