From 93f12e0c74eedd37f21888af03cd6e9eb6537ecb Mon Sep 17 00:00:00 2001 From: Stefano Chiapedi Date: Sat, 6 Aug 2016 11:27:07 +0200 Subject: [PATCH 1/9] improved performance --- pygem/freeform.py | 230 ++++++++++++++++++++++++---------------------- 1 file changed, 119 insertions(+), 111 deletions(-) diff --git a/pygem/freeform.py b/pygem/freeform.py index 81d2655..74527ff 100644 --- a/pygem/freeform.py +++ b/pygem/freeform.py @@ -44,116 +44,125 @@ class FFD(object): - """ - Class that handles the Free Form Deformation on the mesh points. + """ + Class that handles the Free Form Deformation on the mesh points. - :param FFDParameters ffd_parameters: parameters of the Free Form Deformation. - :param numpy.ndarray original_mesh_points: coordinates of the original points of the mesh. + :param FFDParameters ffd_parameters: parameters of the Free Form Deformation. + :param numpy.ndarray original_mesh_points: coordinates of the original points of the mesh. - :cvar FFDParameters parameters: parameters of the Free Form Deformation. - :cvar numpy.ndarray original_mesh_points: coordinates of the original points of the mesh. - The shape is `n_points`-by-3. - :cvar numpy.ndarray modified_mesh_points: coordinates of the points of the deformed mesh. + :cvar FFDParameters parameters: parameters of the Free Form Deformation. + :cvar numpy.ndarray original_mesh_points: coordinates of the original points of the mesh. The shape is `n_points`-by-3. + :cvar numpy.ndarray modified_mesh_points: coordinates of the points of the deformed mesh. + The shape is `n_points`-by-3. - :Example: + :Example: - >>> import pygem.freeform as ffd - >>> import pygem.params as ffdp - >>> import numpy as np + >>> import pygem.freeform as ffd + >>> import pygem.params as ffdp + >>> import numpy as np - >>> ffd_parameters = ffdp.FFDParameters() - >>> ffd_parameters.read_parameters('tests/test_datasets/parameters_test_ffd_sphere.prm') - >>> original_mesh_points = np.load('tests/test_datasets/meshpoints_sphere_orig.npy') - >>> free_form = ffd.FFD(ffd_parameters, original_mesh_points) - >>> free_form.perform() - >>> new_mesh_points = free_form.modified_mesh_points + >>> ffd_parameters = ffdp.FFDParameters() + >>> ffd_parameters.read_parameters('tests/test_datasets/parameters_test_ffd_sphere.prm') + >>> original_mesh_points = np.load('tests/test_datasets/meshpoints_sphere_orig.npy') + >>> free_form = ffd.FFD(ffd_parameters, original_mesh_points) + >>> free_form.perform() + >>> new_mesh_points = free_form.modified_mesh_points """ - def __init__(self, ffd_parameters, original_mesh_points): - self.parameters = ffd_parameters - self.original_mesh_points = original_mesh_points - self.modified_mesh_points = None - - - def perform(self): - """ - This method performs the deformation on the mesh points. After the execution - it sets `self.modified_mesh_points`. - - .. todo:: - In order to improve the performances, we need to perform the FFD only on the points inside - the lattice. - """ - (n_rows_mesh, n_cols_mesh) = self.original_mesh_points.shape - - # translation and then affine transformation - translation = self.parameters.origin_box - - fisical_frame = np.array([self.parameters.position_vertex_1 - translation, \ - self.parameters.position_vertex_2 - translation, \ - self.parameters.position_vertex_3 - translation, \ - [0, 0, 0]]) - reference_frame = np.array([[1, 0, 0], [0, 1, 0], [0, 0, 1], [0, 0, 0]]) - - transformation = at.affine_points_fit(fisical_frame, reference_frame) - inverse_transformation = at.affine_points_fit(reference_frame, fisical_frame) - - # Initialization. In order to exploit the contiguity in memory the following are transposed - (dim_n_mu, dim_m_mu, dim_t_mu) = self.parameters.array_mu_x.shape - bernstein_x = np.zeros((dim_n_mu, n_rows_mesh)) - bernstein_y = np.zeros((dim_m_mu, n_rows_mesh)) - bernstein_z = np.zeros((dim_t_mu, n_rows_mesh)) - shift_mesh_points = np.zeros((n_cols_mesh, n_rows_mesh)) - - reference_frame_mesh_points = self._transform_points(self.original_mesh_points - translation, \ - transformation) - - for i in range(0, dim_n_mu): - aux1 = np.power((1-reference_frame_mesh_points[:, 0]), dim_n_mu-1-i) - aux2 = np.power(reference_frame_mesh_points[:, 0], i) - bernstein_x[i, :] = special.binom(dim_n_mu-1, i) * np.multiply(aux1, aux2) - - for i in range(0, dim_m_mu): - aux1 = np.power((1-reference_frame_mesh_points[:, 1]), dim_m_mu-1-i) - aux2 = np.power(reference_frame_mesh_points[:, 1], i) - bernstein_y[i, :] = special.binom(dim_m_mu-1, i) * np.multiply(aux1, aux2) - - for i in range(0, dim_t_mu): - aux1 = np.power((1-reference_frame_mesh_points[:, 2]), dim_t_mu-1-i) - aux2 = np.power(reference_frame_mesh_points[:, 2], i) - bernstein_z[i, :] = special.binom(dim_t_mu-1, i) * np.multiply(aux1, aux2) - - - aux_x = 0. - aux_y = 0. - aux_z = 0. - for j in range(0, dim_m_mu): - for k in range(0, dim_t_mu): - bernstein_yz = np.multiply(bernstein_y[j, :], bernstein_z[k, :]) - for i in range(0, dim_n_mu): - aux = np.multiply(bernstein_x[i, :], bernstein_yz) - aux_x += aux * self.parameters.array_mu_x[i, j, k] - aux_y += aux * self.parameters.array_mu_y[i, j, k] - aux_z += aux * self.parameters.array_mu_z[i, j, k] - shift_mesh_points[0, :] += aux_x - shift_mesh_points[1, :] += aux_y - shift_mesh_points[2, :] += aux_z - - # Splitting points inside and outside the lattice: TODO not very efficient - for i in range(0, n_rows_mesh): - if (reference_frame_mesh_points[i, 0] < 0) or (reference_frame_mesh_points[i, 1] < 0) \ - or (reference_frame_mesh_points[i, 2] < 0) or (reference_frame_mesh_points[i, 0] > 1) \ - or (reference_frame_mesh_points[i, 1] > 1) or (reference_frame_mesh_points[i, 2] > 1): - shift_mesh_points[:, i] = 0 - - # Here shift_mesh_points needs to be transposed to be summed with reference_frame_mesh_points - self.modified_mesh_points = self._transform_points(np.transpose(shift_mesh_points) + \ - reference_frame_mesh_points, inverse_transformation) + translation - - - @staticmethod - def _transform_points(original_points, transformation): - """ + def __init__(self, ffd_parameters, original_mesh_points): + self.parameters = ffd_parameters + self.original_mesh_points = original_mesh_points + self.modified_mesh_points = None + + + def perform(self): + """ + This method performs the deformation on the mesh points. After the execution + it sets `self.modified_mesh_points`. + + .. todo:: +# In order to improve the performances, we need to perform the FFD only on the points inside +# the lattice. + """ + # translation and then affine transformation + translation = self.parameters.origin_box + + physical_frame = np.array([self.parameters.position_vertex_1 - translation, \ + self.parameters.position_vertex_2 - translation, \ + self.parameters.position_vertex_3 - translation, \ + [0, 0, 0]]) + reference_frame = np.array([[1, 0, 0], [0, 1, 0], [0, 0, 1], [0, 0, 0]]) + + transformation = at.affine_points_fit(physical_frame, reference_frame) + inverse_transformation = at.affine_points_fit(reference_frame, physical_frame) + + box_position_x = sorted([self.parameters.origin_box[0], self.parameters.origin_box[0]+self.parameters.lenght_box_x], key=float) + box_position_y = sorted([self.parameters.origin_box[1], self.parameters.origin_box[1]+self.parameters.lenght_box_y], key=float) + box_position_z = sorted([self.parameters.origin_box[2], self.parameters.origin_box[2]+self.parameters.lenght_box_z], key=float) + mesh_points = self.original_mesh_points[(self.original_mesh_points[:,0] >= box_position_x[0]) & \ + (self.original_mesh_points[:,0] <= box_position_x[1]) & \ + (self.original_mesh_points[:,1] >= box_position_y[0]) & \ + (self.original_mesh_points[:,1] <= box_position_y[1]) & \ + (self.original_mesh_points[:,2] >= box_position_z[0]) & \ + (self.original_mesh_points[:,2] <= box_position_z[1])] + (n_rows_mesh, n_cols_mesh) = mesh_points.shape + + # Initialization. In order to exploit the contiguity in memory the following are transposed + (dim_n_mu, dim_m_mu, dim_t_mu) = self.parameters.array_mu_x.shape + bernstein_x = np.zeros((dim_n_mu, n_rows_mesh)) + bernstein_y = np.zeros((dim_m_mu, n_rows_mesh)) + bernstein_z = np.zeros((dim_t_mu, n_rows_mesh)) + shift_mesh_points = np.zeros((n_cols_mesh, n_rows_mesh)) + + reference_frame_mesh_points = self._transform_points(mesh_points - translation, transformation) + + for i in range(0, dim_n_mu): + aux1 = np.power((1-reference_frame_mesh_points[:, 0]), dim_n_mu-1-i) + aux2 = np.power(reference_frame_mesh_points[:, 0], i) + bernstein_x[i, :] = special.binom(dim_n_mu-1, i) * np.multiply(aux1, aux2) + + for i in range(0, dim_m_mu): + aux1 = np.power((1-reference_frame_mesh_points[:, 1]), dim_m_mu-1-i) + aux2 = np.power(reference_frame_mesh_points[:, 1], i) + bernstein_y[i, :] = special.binom(dim_m_mu-1, i) * np.multiply(aux1, aux2) + + for i in range(0, dim_t_mu): + aux1 = np.power((1-reference_frame_mesh_points[:, 2]), dim_t_mu-1-i) + aux2 = np.power(reference_frame_mesh_points[:, 2], i) + bernstein_z[i, :] = special.binom(dim_t_mu-1, i) * np.multiply(aux1, aux2) + + + aux_x = 0. + aux_y = 0. + aux_z = 0. + for j in range(0, dim_m_mu): + for k in range(0, dim_t_mu): + bernstein_yz = np.multiply(bernstein_y[j, :], bernstein_z[k, :]) + for i in range(0, dim_n_mu): + aux = np.multiply(bernstein_x[i, :], bernstein_yz) + aux_x += aux * self.parameters.array_mu_x[i, j, k] + aux_y += aux * self.parameters.array_mu_y[i, j, k] + aux_z += aux * self.parameters.array_mu_z[i, j, k] + shift_mesh_points[0, :] += aux_x + shift_mesh_points[1, :] += aux_y + shift_mesh_points[2, :] += aux_z + + # Here shift_mesh_points needs to be transposed to be summed with reference_frame_mesh_points + new_mesh_points = self._transform_points(np.transpose(shift_mesh_points) + reference_frame_mesh_points, \ + inverse_transformation) + translation + + self.modified_mesh_points = np.copy(self.original_mesh_points) + self.modified_mesh_points[(self.original_mesh_points[:,0] >= box_position_x[0]) & \ + (self.original_mesh_points[:,0] <= box_position_x[1]) & \ + (self.original_mesh_points[:,1] >= box_position_y[0]) & \ + (self.original_mesh_points[:,1] <= box_position_y[1]) & \ + (self.original_mesh_points[:,2] >= box_position_z[0]) & \ + (self.original_mesh_points[:,2] <= box_position_z[1])] \ + = new_mesh_points + + @staticmethod + def _transform_points(original_points, transformation): + """ This private static method transforms the points according to the affine transformation taken from affine_points_fit method. @@ -162,11 +171,10 @@ def _transform_points(original_points, transformation): :return: modified_points: coordinates of the modified points. :rtype: numpy.ndarray - """ - n_rows_mesh, n_cols_mesh = original_points.shape - modified_points = np.zeros((n_rows_mesh, n_cols_mesh)) - - for i in range(0, n_rows_mesh): - modified_points[i, :] = transformation(original_points[i]) - - return modified_points + """ + n_rows_mesh, n_cols_mesh = original_points.shape + modified_points = np.zeros((n_rows_mesh, n_cols_mesh)) + for i in range(0, n_rows_mesh): + modified_points[i, :] = transformation(original_points[i]) + + return modified_points From 6895ee4746c63cebb1b6b2802ae973c683f362e0 Mon Sep 17 00:00:00 2001 From: Stefano Chiapedi Date: Sat, 6 Aug 2016 11:27:50 +0200 Subject: [PATCH 2/9] fixed bug and improved flexibility --- pygem/unvhandler.py | 91 ++++++++++++++++++++------------------------- 1 file changed, 41 insertions(+), 50 deletions(-) diff --git a/pygem/unvhandler.py b/pygem/unvhandler.py index 8262b4b..9f21a63 100644 --- a/pygem/unvhandler.py +++ b/pygem/unvhandler.py @@ -34,43 +34,28 @@ def parse(self, filename): self.infile = filename - with open(self.infile, 'r') as input_file: - nline = 0 - while True: - line = input_file.readline() - nline += 1 - if not line: - break - if line.startswith(' -1'): - section_id = input_file.readline().strip() - nline += 1 - if section_id == '2411': - count = 0 - while not input_file.readline().startswith(' -1'): - count += 1 - start_line = nline + 2 - last_line = start_line + count - else: - while not input_file.readline().startswith(' -1'): - nline += 1 - - n_points = count/2 - mesh_points = np.zeros(shape=(n_points, 3)) - - nline = 0 - i = 0 - with open(self.infile, 'r') as input_file: - for line in input_file: - nline += 1 - if nline % 2 == 1 and start_line < nline < last_line: - line = line.strip() - j = 0 - for number in line.split(): - mesh_points[i][j] = float(number) - j += 1 - i += 1 - - return mesh_points + index = -9 + mesh_points = [] + with open(self.infile, 'r') as input_file: + for num, line in enumerate(input_file): + if line.startswith(' 2411'): + index = num + if num == index + 2: + if line.startswith(' -1'): + break + else: + line = line.replace('D', 'E') + l = [] + for t in line.split(): + try: + l.append(float(t)) + except ValueError: + pass + mesh_points.append(l) + index = num + mesh_points = np.array(mesh_points) + + return mesh_points def write(self, mesh_points, filename): @@ -89,16 +74,22 @@ def write(self, mesh_points, filename): self.outfile = filename - n_points = mesh_points.shape[0] - nrow = 0 - i = 0 - with open(self.infile, 'r') as input_file, open(self.outfile, 'w') as output_file: - for line in input_file: - nrow += 1 - if nrow % 2 == 1 and 20 < nrow <= (20 + n_points * 2): - for j in range(0, 3): - output_file.write(' ' + str(mesh_points[i][j])) - output_file.write('\n') - i += 1 - elif nrow > 17: - output_file.write(line) + index = -9 + i = 0 + with open(self.outfile, 'w') as output_file: + with open(self.infile, 'r') as input_file: + for num, line in enumerate(input_file): + if line.startswith(' 2411'): + index = num + if num == index + 2: + if line.startswith(' -1'): + index = -9 + output_file.write(line) + else: + for j in range(0, 3): + output_file.write(3*' ' + '{:.16E}'.format(mesh_points[i][j])) + output_file.write('\n') + i += 1 + index = num + else: + output_file.write(line) From 2a41397b55d25dd8575ddd0701e9b52ca9840042 Mon Sep 17 00:00:00 2001 From: Stefano Chiapedi Date: Sun, 7 Aug 2016 00:18:46 +0200 Subject: [PATCH 3/9] converted indentations with tabs --- pygem/freeform.py | 246 +++++++++++++++++++++++----------------------- 1 file changed, 123 insertions(+), 123 deletions(-) diff --git a/pygem/freeform.py b/pygem/freeform.py index 74527ff..969b178 100644 --- a/pygem/freeform.py +++ b/pygem/freeform.py @@ -44,125 +44,124 @@ class FFD(object): - """ - Class that handles the Free Form Deformation on the mesh points. + """ + Class that handles the Free Form Deformation on the mesh points. - :param FFDParameters ffd_parameters: parameters of the Free Form Deformation. - :param numpy.ndarray original_mesh_points: coordinates of the original points of the mesh. + :param FFDParameters ffd_parameters: parameters of the Free Form Deformation. + :param numpy.ndarray original_mesh_points: coordinates of the original points of the mesh. - :cvar FFDParameters parameters: parameters of the Free Form Deformation. - :cvar numpy.ndarray original_mesh_points: coordinates of the original points of the mesh. + :cvar FFDParameters parameters: parameters of the Free Form Deformation. + :cvar numpy.ndarray original_mesh_points: coordinates of the original points of the mesh. The shape is `n_points`-by-3. - :cvar numpy.ndarray modified_mesh_points: coordinates of the points of the deformed mesh. - The shape is `n_points`-by-3. - - :Example: - - >>> import pygem.freeform as ffd - >>> import pygem.params as ffdp - >>> import numpy as np - - >>> ffd_parameters = ffdp.FFDParameters() - >>> ffd_parameters.read_parameters('tests/test_datasets/parameters_test_ffd_sphere.prm') - >>> original_mesh_points = np.load('tests/test_datasets/meshpoints_sphere_orig.npy') - >>> free_form = ffd.FFD(ffd_parameters, original_mesh_points) - >>> free_form.perform() - >>> new_mesh_points = free_form.modified_mesh_points - """ - def __init__(self, ffd_parameters, original_mesh_points): - self.parameters = ffd_parameters - self.original_mesh_points = original_mesh_points - self.modified_mesh_points = None - - - def perform(self): - """ - This method performs the deformation on the mesh points. After the execution - it sets `self.modified_mesh_points`. - - .. todo:: -# In order to improve the performances, we need to perform the FFD only on the points inside -# the lattice. - """ - # translation and then affine transformation - translation = self.parameters.origin_box - - physical_frame = np.array([self.parameters.position_vertex_1 - translation, \ - self.parameters.position_vertex_2 - translation, \ - self.parameters.position_vertex_3 - translation, \ - [0, 0, 0]]) - reference_frame = np.array([[1, 0, 0], [0, 1, 0], [0, 0, 1], [0, 0, 0]]) - - transformation = at.affine_points_fit(physical_frame, reference_frame) - inverse_transformation = at.affine_points_fit(reference_frame, physical_frame) - - box_position_x = sorted([self.parameters.origin_box[0], self.parameters.origin_box[0]+self.parameters.lenght_box_x], key=float) - box_position_y = sorted([self.parameters.origin_box[1], self.parameters.origin_box[1]+self.parameters.lenght_box_y], key=float) - box_position_z = sorted([self.parameters.origin_box[2], self.parameters.origin_box[2]+self.parameters.lenght_box_z], key=float) - mesh_points = self.original_mesh_points[(self.original_mesh_points[:,0] >= box_position_x[0]) & \ - (self.original_mesh_points[:,0] <= box_position_x[1]) & \ - (self.original_mesh_points[:,1] >= box_position_y[0]) & \ - (self.original_mesh_points[:,1] <= box_position_y[1]) & \ - (self.original_mesh_points[:,2] >= box_position_z[0]) & \ - (self.original_mesh_points[:,2] <= box_position_z[1])] - (n_rows_mesh, n_cols_mesh) = mesh_points.shape - - # Initialization. In order to exploit the contiguity in memory the following are transposed - (dim_n_mu, dim_m_mu, dim_t_mu) = self.parameters.array_mu_x.shape - bernstein_x = np.zeros((dim_n_mu, n_rows_mesh)) - bernstein_y = np.zeros((dim_m_mu, n_rows_mesh)) - bernstein_z = np.zeros((dim_t_mu, n_rows_mesh)) - shift_mesh_points = np.zeros((n_cols_mesh, n_rows_mesh)) - - reference_frame_mesh_points = self._transform_points(mesh_points - translation, transformation) - - for i in range(0, dim_n_mu): - aux1 = np.power((1-reference_frame_mesh_points[:, 0]), dim_n_mu-1-i) - aux2 = np.power(reference_frame_mesh_points[:, 0], i) - bernstein_x[i, :] = special.binom(dim_n_mu-1, i) * np.multiply(aux1, aux2) - - for i in range(0, dim_m_mu): - aux1 = np.power((1-reference_frame_mesh_points[:, 1]), dim_m_mu-1-i) - aux2 = np.power(reference_frame_mesh_points[:, 1], i) - bernstein_y[i, :] = special.binom(dim_m_mu-1, i) * np.multiply(aux1, aux2) - - for i in range(0, dim_t_mu): - aux1 = np.power((1-reference_frame_mesh_points[:, 2]), dim_t_mu-1-i) - aux2 = np.power(reference_frame_mesh_points[:, 2], i) - bernstein_z[i, :] = special.binom(dim_t_mu-1, i) * np.multiply(aux1, aux2) - - - aux_x = 0. - aux_y = 0. - aux_z = 0. - for j in range(0, dim_m_mu): - for k in range(0, dim_t_mu): - bernstein_yz = np.multiply(bernstein_y[j, :], bernstein_z[k, :]) - for i in range(0, dim_n_mu): - aux = np.multiply(bernstein_x[i, :], bernstein_yz) - aux_x += aux * self.parameters.array_mu_x[i, j, k] - aux_y += aux * self.parameters.array_mu_y[i, j, k] - aux_z += aux * self.parameters.array_mu_z[i, j, k] - shift_mesh_points[0, :] += aux_x - shift_mesh_points[1, :] += aux_y - shift_mesh_points[2, :] += aux_z - - # Here shift_mesh_points needs to be transposed to be summed with reference_frame_mesh_points - new_mesh_points = self._transform_points(np.transpose(shift_mesh_points) + reference_frame_mesh_points, \ - inverse_transformation) + translation - - self.modified_mesh_points = np.copy(self.original_mesh_points) - self.modified_mesh_points[(self.original_mesh_points[:,0] >= box_position_x[0]) & \ - (self.original_mesh_points[:,0] <= box_position_x[1]) & \ - (self.original_mesh_points[:,1] >= box_position_y[0]) & \ - (self.original_mesh_points[:,1] <= box_position_y[1]) & \ - (self.original_mesh_points[:,2] >= box_position_z[0]) & \ - (self.original_mesh_points[:,2] <= box_position_z[1])] \ - = new_mesh_points - - @staticmethod - def _transform_points(original_points, transformation): - """ + :cvar numpy.ndarray modified_mesh_points: coordinates of the points of the deformed mesh. + The shape is `n_points`-by-3. + + :Example: + + >>> import pygem.freeform as ffd + >>> import pygem.params as ffdp + >>> import numpy as np + + >>> ffd_parameters = ffdp.FFDParameters() + >>> ffd_parameters.read_parameters('tests/test_datasets/parameters_test_ffd_sphere.prm') + >>> original_mesh_points = np.load('tests/test_datasets/meshpoints_sphere_orig.npy') + >>> free_form = ffd.FFD(ffd_parameters, original_mesh_points) + >>> free_form.perform() + >>> new_mesh_points = free_form.modified_mesh_points + """ + def __init__(self, ffd_parameters, original_mesh_points): + self.parameters = ffd_parameters + self.original_mesh_points = original_mesh_points + self.modified_mesh_points = None + + + def perform(self): + """ + This method performs the deformation on the mesh points. After the execution + it sets `self.modified_mesh_points`. + + .. todo:: + + """ + # translation and then affine transformation + translation = self.parameters.origin_box + + physical_frame = np.array([self.parameters.position_vertex_1 - translation, \ + self.parameters.position_vertex_2 - translation, \ + self.parameters.position_vertex_3 - translation, \ + [0, 0, 0]]) + reference_frame = np.array([[1, 0, 0], [0, 1, 0], [0, 0, 1], [0, 0, 0]]) + + transformation = at.affine_points_fit(physical_frame, reference_frame) + inverse_transformation = at.affine_points_fit(reference_frame, physical_frame) + + box_position_x = sorted([self.parameters.origin_box[0], self.parameters.origin_box[0]+self.parameters.lenght_box_x], key=float) + box_position_y = sorted([self.parameters.origin_box[1], self.parameters.origin_box[1]+self.parameters.lenght_box_y], key=float) + box_position_z = sorted([self.parameters.origin_box[2], self.parameters.origin_box[2]+self.parameters.lenght_box_z], key=float) + mesh_points = self.original_mesh_points[(self.original_mesh_points[:,0] >= box_position_x[0]) & \ + (self.original_mesh_points[:,0] <= box_position_x[1]) & \ + (self.original_mesh_points[:,1] >= box_position_y[0]) & \ + (self.original_mesh_points[:,1] <= box_position_y[1]) & \ + (self.original_mesh_points[:,2] >= box_position_z[0]) & \ + (self.original_mesh_points[:,2] <= box_position_z[1])] + (n_rows_mesh, n_cols_mesh) = mesh_points.shape + + # Initialization. In order to exploit the contiguity in memory the following are transposed + (dim_n_mu, dim_m_mu, dim_t_mu) = self.parameters.array_mu_x.shape + bernstein_x = np.zeros((dim_n_mu, n_rows_mesh)) + bernstein_y = np.zeros((dim_m_mu, n_rows_mesh)) + bernstein_z = np.zeros((dim_t_mu, n_rows_mesh)) + shift_mesh_points = np.zeros((n_cols_mesh, n_rows_mesh)) + + reference_frame_mesh_points = self._transform_points(mesh_points - translation, transformation) + + for i in range(0, dim_n_mu): + aux1 = np.power((1-reference_frame_mesh_points[:, 0]), dim_n_mu-1-i) + aux2 = np.power(reference_frame_mesh_points[:, 0], i) + bernstein_x[i, :] = special.binom(dim_n_mu-1, i) * np.multiply(aux1, aux2) + + for i in range(0, dim_m_mu): + aux1 = np.power((1-reference_frame_mesh_points[:, 1]), dim_m_mu-1-i) + aux2 = np.power(reference_frame_mesh_points[:, 1], i) + bernstein_y[i, :] = special.binom(dim_m_mu-1, i) * np.multiply(aux1, aux2) + + for i in range(0, dim_t_mu): + aux1 = np.power((1-reference_frame_mesh_points[:, 2]), dim_t_mu-1-i) + aux2 = np.power(reference_frame_mesh_points[:, 2], i) + bernstein_z[i, :] = special.binom(dim_t_mu-1, i) * np.multiply(aux1, aux2) + + + aux_x = 0. + aux_y = 0. + aux_z = 0. + for j in range(0, dim_m_mu): + for k in range(0, dim_t_mu): + bernstein_yz = np.multiply(bernstein_y[j, :], bernstein_z[k, :]) + for i in range(0, dim_n_mu): + aux = np.multiply(bernstein_x[i, :], bernstein_yz) + aux_x += aux * self.parameters.array_mu_x[i, j, k] + aux_y += aux * self.parameters.array_mu_y[i, j, k] + aux_z += aux * self.parameters.array_mu_z[i, j, k] + shift_mesh_points[0, :] += aux_x + shift_mesh_points[1, :] += aux_y + shift_mesh_points[2, :] += aux_z + + # Here shift_mesh_points needs to be transposed to be summed with reference_frame_mesh_points + new_mesh_points = self._transform_points(np.transpose(shift_mesh_points) + reference_frame_mesh_points, \ + inverse_transformation) + translation + + self.modified_mesh_points = np.copy(self.original_mesh_points) + self.modified_mesh_points[(self.original_mesh_points[:,0] >= box_position_x[0]) & \ + (self.original_mesh_points[:,0] <= box_position_x[1]) & \ + (self.original_mesh_points[:,1] >= box_position_y[0]) & \ + (self.original_mesh_points[:,1] <= box_position_y[1]) & \ + (self.original_mesh_points[:,2] >= box_position_z[0]) & \ + (self.original_mesh_points[:,2] <= box_position_z[1])] \ + = new_mesh_points + + @staticmethod + def _transform_points(original_points, transformation): + """ This private static method transforms the points according to the affine transformation taken from affine_points_fit method. @@ -171,10 +170,11 @@ def _transform_points(original_points, transformation): :return: modified_points: coordinates of the modified points. :rtype: numpy.ndarray - """ - n_rows_mesh, n_cols_mesh = original_points.shape - modified_points = np.zeros((n_rows_mesh, n_cols_mesh)) - for i in range(0, n_rows_mesh): - modified_points[i, :] = transformation(original_points[i]) - - return modified_points + """ + n_rows_mesh, n_cols_mesh = original_points.shape + modified_points = np.zeros((n_rows_mesh, n_cols_mesh)) + + for i in range(0, n_rows_mesh): + modified_points[i, :] = transformation(original_points[i]) + + return modified_points From 69f7694965a4f35a79f935fb4c1743b2d651d936 Mon Sep 17 00:00:00 2001 From: Stefano Chiapedi Date: Sun, 7 Aug 2016 00:19:30 +0200 Subject: [PATCH 4/9] converted indentations with tabs --- pygem/unvhandler.py | 80 ++++++++++++++++++++++----------------------- 1 file changed, 40 insertions(+), 40 deletions(-) diff --git a/pygem/unvhandler.py b/pygem/unvhandler.py index 9f21a63..9c35ad4 100644 --- a/pygem/unvhandler.py +++ b/pygem/unvhandler.py @@ -34,28 +34,28 @@ def parse(self, filename): self.infile = filename - index = -9 - mesh_points = [] - with open(self.infile, 'r') as input_file: - for num, line in enumerate(input_file): - if line.startswith(' 2411'): - index = num - if num == index + 2: - if line.startswith(' -1'): - break - else: - line = line.replace('D', 'E') - l = [] - for t in line.split(): - try: - l.append(float(t)) - except ValueError: - pass - mesh_points.append(l) - index = num - mesh_points = np.array(mesh_points) + index = -9 + mesh_points = [] + with open(self.infile, 'r') as input_file: + for num, line in enumerate(input_file): + if line.startswith(' 2411'): + index = num + if num == index + 2: + if line.startswith(' -1'): + break + else: + line = line.replace('D', 'E') + l = [] + for t in line.split(): + try: + l.append(float(t)) + except ValueError: + pass + mesh_points.append(l) + index = num + mesh_points = np.array(mesh_points) - return mesh_points + return mesh_points def write(self, mesh_points, filename): @@ -74,22 +74,22 @@ def write(self, mesh_points, filename): self.outfile = filename - index = -9 - i = 0 - with open(self.outfile, 'w') as output_file: - with open(self.infile, 'r') as input_file: - for num, line in enumerate(input_file): - if line.startswith(' 2411'): - index = num - if num == index + 2: - if line.startswith(' -1'): - index = -9 - output_file.write(line) - else: - for j in range(0, 3): - output_file.write(3*' ' + '{:.16E}'.format(mesh_points[i][j])) - output_file.write('\n') - i += 1 - index = num - else: - output_file.write(line) + index = -9 + i = 0 + with open(self.outfile, 'w') as output_file: + with open(self.infile, 'r') as input_file: + for num, line in enumerate(input_file): + if line.startswith(' 2411'): + index = num + if num == index + 2: + if line.startswith(' -1'): + index = -9 + output_file.write(line) + else: + for j in range(0, 3): + output_file.write(3*' ' + '{:.16E}'.format(mesh_points[i][j])) + output_file.write('\n') + i += 1 + index = num + else: + output_file.write(line) From 2466494b4b3a29e5fc03814fc1603fe11198afb3 Mon Sep 17 00:00:00 2001 From: Stefano Chiapedi Date: Sun, 7 Aug 2016 18:48:37 +0200 Subject: [PATCH 5/9] fixed error with box rotation --- pygem/freeform.py | 68 ++++++++++++++++++++++++----------------------- 1 file changed, 35 insertions(+), 33 deletions(-) diff --git a/pygem/freeform.py b/pygem/freeform.py index 969b178..f991635 100644 --- a/pygem/freeform.py +++ b/pygem/freeform.py @@ -86,24 +86,25 @@ def perform(self): # translation and then affine transformation translation = self.parameters.origin_box - physical_frame = np.array([self.parameters.position_vertex_1 - translation, \ - self.parameters.position_vertex_2 - translation, \ - self.parameters.position_vertex_3 - translation, \ + physical_frame = np.array([self.parameters.position_vertex_1 - translation, + self.parameters.position_vertex_2 - translation, + self.parameters.position_vertex_3 - translation, [0, 0, 0]]) reference_frame = np.array([[1, 0, 0], [0, 1, 0], [0, 0, 1], [0, 0, 0]]) transformation = at.affine_points_fit(physical_frame, reference_frame) inverse_transformation = at.affine_points_fit(reference_frame, physical_frame) - - box_position_x = sorted([self.parameters.origin_box[0], self.parameters.origin_box[0]+self.parameters.lenght_box_x], key=float) - box_position_y = sorted([self.parameters.origin_box[1], self.parameters.origin_box[1]+self.parameters.lenght_box_y], key=float) - box_position_z = sorted([self.parameters.origin_box[2], self.parameters.origin_box[2]+self.parameters.lenght_box_z], key=float) - mesh_points = self.original_mesh_points[(self.original_mesh_points[:,0] >= box_position_x[0]) & \ - (self.original_mesh_points[:,0] <= box_position_x[1]) & \ - (self.original_mesh_points[:,1] >= box_position_y[0]) & \ - (self.original_mesh_points[:,1] <= box_position_y[1]) & \ - (self.original_mesh_points[:,2] >= box_position_z[0]) & \ - (self.original_mesh_points[:,2] <= box_position_z[1])] + + # apply transformation to original mesh points + reference_frame_mesh_points = self._transform_points(self.original_mesh_points - translation, transformation) + + # select mesh points inside bounding box + mesh_points = reference_frame_mesh_points[(reference_frame_mesh_points[:,0] >= 0.) & + (reference_frame_mesh_points[:,0] <= 1.) & + (reference_frame_mesh_points[:,1] >= 0.) & + (reference_frame_mesh_points[:,1] <= 1.) & + (reference_frame_mesh_points[:,2] >= 0.) & + (reference_frame_mesh_points[:,2] <= 1.)] (n_rows_mesh, n_cols_mesh) = mesh_points.shape # Initialization. In order to exploit the contiguity in memory the following are transposed @@ -113,21 +114,19 @@ def perform(self): bernstein_z = np.zeros((dim_t_mu, n_rows_mesh)) shift_mesh_points = np.zeros((n_cols_mesh, n_rows_mesh)) - reference_frame_mesh_points = self._transform_points(mesh_points - translation, transformation) - for i in range(0, dim_n_mu): - aux1 = np.power((1-reference_frame_mesh_points[:, 0]), dim_n_mu-1-i) - aux2 = np.power(reference_frame_mesh_points[:, 0], i) + aux1 = np.power((1-mesh_points[:, 0]), dim_n_mu-1-i) + aux2 = np.power(mesh_points[:, 0], i) bernstein_x[i, :] = special.binom(dim_n_mu-1, i) * np.multiply(aux1, aux2) for i in range(0, dim_m_mu): - aux1 = np.power((1-reference_frame_mesh_points[:, 1]), dim_m_mu-1-i) - aux2 = np.power(reference_frame_mesh_points[:, 1], i) + aux1 = np.power((1-mesh_points[:, 1]), dim_m_mu-1-i) + aux2 = np.power(mesh_points[:, 1], i) bernstein_y[i, :] = special.binom(dim_m_mu-1, i) * np.multiply(aux1, aux2) for i in range(0, dim_t_mu): - aux1 = np.power((1-reference_frame_mesh_points[:, 2]), dim_t_mu-1-i) - aux2 = np.power(reference_frame_mesh_points[:, 2], i) + aux1 = np.power((1-mesh_points[:, 2]), dim_t_mu-1-i) + aux2 = np.power(mesh_points[:, 2], i) bernstein_z[i, :] = special.binom(dim_t_mu-1, i) * np.multiply(aux1, aux2) @@ -146,18 +145,21 @@ def perform(self): shift_mesh_points[1, :] += aux_y shift_mesh_points[2, :] += aux_z - # Here shift_mesh_points needs to be transposed to be summed with reference_frame_mesh_points - new_mesh_points = self._transform_points(np.transpose(shift_mesh_points) + reference_frame_mesh_points, \ - inverse_transformation) + translation - - self.modified_mesh_points = np.copy(self.original_mesh_points) - self.modified_mesh_points[(self.original_mesh_points[:,0] >= box_position_x[0]) & \ - (self.original_mesh_points[:,0] <= box_position_x[1]) & \ - (self.original_mesh_points[:,1] >= box_position_y[0]) & \ - (self.original_mesh_points[:,1] <= box_position_y[1]) & \ - (self.original_mesh_points[:,2] >= box_position_z[0]) & \ - (self.original_mesh_points[:,2] <= box_position_z[1])] \ - = new_mesh_points + # shift_mesh_points needs to be transposed to be summed with mesh_points + # apply inverse transformation to shifted mesh points + new_mesh_points = self._transform_points(np.transpose(shift_mesh_points) + + mesh_points, inverse_transformation) + \ + translation + + # merge non-shifted mesh points with shifted ones + self.modified_mesh_points = np.copy(self.original_mesh_points) + self.modified_mesh_points[(reference_frame_mesh_points[:,0] >= 0.) & + (reference_frame_mesh_points[:,0] <= 1.) & + (reference_frame_mesh_points[:,1] >= 0.) & + (reference_frame_mesh_points[:,1] <= 1.) & + (reference_frame_mesh_points[:,2] >= 0.) & + (reference_frame_mesh_points[:,2] <= 1.)] \ + = new_mesh_points @staticmethod def _transform_points(original_points, transformation): From 1443987f813897139a868edd5c7fb0b1474fbc28 Mon Sep 17 00:00:00 2001 From: Stefano Chiapedi Date: Sun, 7 Aug 2016 18:49:37 +0200 Subject: [PATCH 6/9] reconverted tabs in spaces in startswith command --- pygem/unvhandler.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/pygem/unvhandler.py b/pygem/unvhandler.py index 9c35ad4..cd1773c 100644 --- a/pygem/unvhandler.py +++ b/pygem/unvhandler.py @@ -41,7 +41,7 @@ def parse(self, filename): if line.startswith(' 2411'): index = num if num == index + 2: - if line.startswith(' -1'): + if line.startswith(' -1'): break else: line = line.replace('D', 'E') @@ -82,7 +82,7 @@ def write(self, mesh_points, filename): if line.startswith(' 2411'): index = num if num == index + 2: - if line.startswith(' -1'): + if line.startswith(' -1'): index = -9 output_file.write(line) else: From fd7fe28c34da2283a50c1d7600fe19c1735c011c Mon Sep 17 00:00:00 2001 From: Stefano Chiapedi Date: Sun, 7 Aug 2016 19:51:03 +0200 Subject: [PATCH 7/9] new write comparison --- tests/test_unvhandler.py | 8 +------- 1 file changed, 1 insertion(+), 7 deletions(-) diff --git a/tests/test_unvhandler.py b/tests/test_unvhandler.py index 42f6b27..e6e99a3 100644 --- a/tests/test_unvhandler.py +++ b/tests/test_unvhandler.py @@ -116,15 +116,9 @@ def test_unv_write_outfile(self): def test_unv_write_comparison(self): unv_handler = uh.UnvHandler() mesh_points = unv_handler.parse('tests/test_datasets/test_square.unv') - mesh_points[0][0] = 2.2 - mesh_points[5][1] = 4.3 - mesh_points[9][2] = 0.5 - mesh_points[45][0] = 7.2 - mesh_points[132][1] = -1.2 - mesh_points[255][2] = -3.6 outfilename = 'tests/test_datasets/test_square_out.unv' - outfilename_expected = 'tests/test_datasets/test_square_out_true.unv' + outfilename_expected = 'tests/test_datasets/test_square.unv' unv_handler.write(mesh_points, outfilename) self.assertTrue(filecmp.cmp(outfilename, outfilename_expected)) From ceed71b38558737ca0025855ea150a2135bad42c Mon Sep 17 00:00:00 2001 From: Stefano Chiapedi Date: Sun, 7 Aug 2016 20:04:48 +0200 Subject: [PATCH 8/9] test_unvhandler changed --- pygem/unvhandler.py | 1 + 1 file changed, 1 insertion(+) diff --git a/pygem/unvhandler.py b/pygem/unvhandler.py index cd1773c..0f9af98 100644 --- a/pygem/unvhandler.py +++ b/pygem/unvhandler.py @@ -93,3 +93,4 @@ def write(self, mesh_points, filename): index = num else: output_file.write(line) + From f61e82c86a54d55c472f4d5b4aa750a7ecc2f0fc Mon Sep 17 00:00:00 2001 From: Stefano Chiapedi Date: Sun, 7 Aug 2016 20:07:34 +0200 Subject: [PATCH 9/9] test_unvhandler changed --- pygem/freeform.py | 1 + 1 file changed, 1 insertion(+) diff --git a/pygem/freeform.py b/pygem/freeform.py index f991635..45176fe 100644 --- a/pygem/freeform.py +++ b/pygem/freeform.py @@ -180,3 +180,4 @@ def _transform_points(original_points, transformation): modified_points[i, :] = transformation(original_points[i]) return modified_points +