From c611ff08549ebe51a1af94840f4f809203440624 Mon Sep 17 00:00:00 2001 From: Bernardo Ferreira Date: Thu, 26 Oct 2023 16:22:49 -0400 Subject: [PATCH] PRv1.0.4 (#11) * Updated version. * Fixed bug when building first Piola-Kirchhoff stress tensor for output. * Fixed execution error under finite strains. - Avoid state update kinematic post-processing when state update fails; * Reduced overhead cost of matrix condensation. * Fixed bug when writting voxels output file. * Fixed bug when checking adaptivity input parameters. * Fixed clusters strains initial iterative guess under clustering adaptivity. --- VERSION | 2 +- docs/source/conf.py | 2 +- src/cratepy/ioput/info.py | 2 +- .../ioput/miscoutputfiles/voxelsoutput.py | 13 +-- src/cratepy/ioput/readprocedures.py | 4 +- src/cratepy/material/materialmodeling.py | 4 +- src/cratepy/online/crom/asca.py | 79 ++++++++++++++++--- src/cratepy/tensor/matrixoperations.py | 44 ++++++----- 8 files changed, 108 insertions(+), 42 deletions(-) diff --git a/VERSION b/VERSION index 21e8796a..ee90284c 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -1.0.3 +1.0.4 diff --git a/docs/source/conf.py b/docs/source/conf.py index 74050cd1..5fc4a885 100644 --- a/docs/source/conf.py +++ b/docs/source/conf.py @@ -23,7 +23,7 @@ author = 'Bernardo Ferreira' copyright = '2020, Bernardo Ferreira' version = '1.0' -release = '1.0.3' +release = '1.0.4' # -- General configuration ---------------------------------------------------- diff --git a/src/cratepy/ioput/info.py b/src/cratepy/ioput/info.py index 1582bd3c..e4fd97c0 100644 --- a/src/cratepy/ioput/info.py +++ b/src/cratepy/ioput/info.py @@ -93,7 +93,7 @@ def displayinfo(code, *args, **kwargs): elif code == '0': arguments = \ ['CRATE - Clustering-based Nonlinear Analysis of Materials', - 'Created by Bernardo P. Ferreira', 'Release 1.0.3 (Jul 2023)'] \ + 'Created by Bernardo P. Ferreira', 'Release 1.0.4 (Oct 2023)'] \ + 2*[args[0], ] + list(args[1:3]) info = tuple(arguments) template = '\n' + colorama.Fore.WHITE + tilde_line \ diff --git a/src/cratepy/ioput/miscoutputfiles/voxelsoutput.py b/src/cratepy/ioput/miscoutputfiles/voxelsoutput.py index 5930d888..845db1b2 100644 --- a/src/cratepy/ioput/miscoutputfiles/voxelsoutput.py +++ b/src/cratepy/ioput/miscoutputfiles/voxelsoutput.py @@ -46,7 +46,8 @@ class VoxelsOutput: ------- init_voxels_output_file(self, crve) Open output file and write file header. - write_voxels_output_file(self, n_dim, comp_order, crve, clusters_state) + write_voxels_output_file(self, n_dim, comp_order, crve, clusters_state, \ + clusters_def_gradient_mf) Write output file. rewind_file(self, rewind_inc) Rewind output file. @@ -103,7 +104,7 @@ def init_voxels_output_file(self, crve): open(self._voxout_file_path, 'w').writelines(write_list) # ------------------------------------------------------------------------- def write_voxels_output_file(self, n_dim, comp_order, crve, - clusters_state): + clusters_state, clusters_def_gradient_mf): """Write output file. Parameters @@ -117,6 +118,9 @@ def write_voxels_output_file(self, n_dim, comp_order, crve, clusters_state : dict Material constitutive model state variables (item, dict) associated to each material cluster (key, str). + clusters_def_gradient_mf : dict + Deformation gradient (item, numpy.ndarray (1d)) associated with + each material cluster (key, str), stored in matricial form. """ # ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ # Instantiate factory of voxels arrays @@ -138,8 +142,7 @@ def write_voxels_output_file(self, n_dim, comp_order, crve, outvar_n_dims = self._output_variables[i][1] # Get material-related output variable voxels array list array_vox_list = voxels_array_factory.build_voxels_array( - crve, outvar, clusters_state) - + crve, outvar, clusters_state, clusters_def_gradient_mf) # Loop over material-related output variable dimensions for var_dim in range(outvar_n_dims): # Store material-related output variable dimension in voxels @@ -298,7 +301,7 @@ def build_voxels_array(self, crve, csbvar, clusters_state, # Build first Piola-Kirchhoff stress tensor first_piola_stress = mop.get_tensor_from_mf( first_piola_stress_mf, self._n_dim, - self._comp_order_sym) + self._comp_order_nsym) # Compute Cauchy stress tensor cauchy_stress = cauchy_from_first_piola( def_gradient, first_piola_stress) diff --git a/src/cratepy/ioput/readprocedures.py b/src/cratepy/ioput/readprocedures.py index f0cdfdc4..d7ae09bc 100644 --- a/src/cratepy/ioput/readprocedures.py +++ b/src/cratepy/ioput/readprocedures.py @@ -1649,7 +1649,7 @@ def read_cluster_analysis_scheme(file, file_path, keyword, material_phases, # Get parameter value value = get_formatted_parameter(parameter, line[1]) # Store adaptivity criterion parameter - if isinstance(value, oacp[parameter]): + if isinstance(value, type(oacp[parameter])): adapt_criterion_data[mat_phase][parameter] = \ get_formatted_parameter(parameter, type(oacp[parameter])(line[1])) @@ -1665,7 +1665,7 @@ def read_cluster_analysis_scheme(file, file_path, keyword, material_phases, # Get parameter value value = get_formatted_parameter(parameter, line[1]) # Store adaptivity type parameter - if isinstance(value, oatp[parameter]): + if isinstance(value, type(oatp[parameter])): adaptivity_type[mat_phase][parameter] = \ get_formatted_parameter(parameter, type(oatp[parameter])(line[1])) diff --git a/src/cratepy/material/materialmodeling.py b/src/cratepy/material/materialmodeling.py index 99595855..a2a339fa 100644 --- a/src/cratepy/material/materialmodeling.py +++ b/src/cratepy/material/materialmodeling.py @@ -1023,7 +1023,9 @@ def _material_su_interface(strain_formulation, problem_type, constitutive_model.state_update(inc_strain, state_variables_old) # ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ # Compute Cauchy stress tensor and material consistent tangent modulus - if strain_formulation == 'finite' and strain_type == 'finite-kinext': + if (not state_variables['is_su_fail'] + and strain_formulation == 'finite' + and strain_type == 'finite-kinext'): # Get Kirchhoff stress tensor (matricial form) kirchhoff_stress_mf = state_variables['stress_mf'] # Build Kirchhoff stress tensor diff --git a/src/cratepy/online/crom/asca.py b/src/cratepy/online/crom/asca.py index b6149526..59b289a0 100644 --- a/src/cratepy/online/crom/asca.py +++ b/src/cratepy/online/crom/asca.py @@ -135,7 +135,7 @@ class ASCA: is_vtk_output=False, \ vtk_data=None, is_voxels_output=False) Solve clustering-based reduced-order equilibrium problem. - _init_global_strain_mf(self, mode='last_converged') + _init_global_strain_mf(self, crve, material_state, mode='last_converged') Set clusters strains initial iterative guess. _init_global_inc_strain_mf(self, n_total_clusters, mode='last_converged') Set clusters incremental strains initial iterative guess. @@ -399,16 +399,16 @@ def solve_equilibrium_problem(self, crve, material_state, mac_load, is_adapt_repeat_inc = True # ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ # Set clustering adaptivity frequency default - if self._clust_adapt_freq is None: - self._clust_adapt_freq = {mat_phase: 1 for mat_phase - in crve.get_adapt_material_phases()} + if clust_adapt_freq is None: + clust_adapt_freq = {mat_phase: 1 for mat_phase + in crve.get_adapt_material_phases()} # ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ # Initialize clustering adaptivity manager adaptivity_manager = AdaptivityManager( self._strain_formulation, self._problem_type, crve.get_adapt_material_phases(), crve.get_phase_clusters(), crve.get_adaptivity_control_feature(), - crve.get_adapt_criterion_data(), self._clust_adapt_freq) + crve.get_adapt_criterion_data(), clust_adapt_freq) # ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ # Initialize increment rewinder manager if is_solution_rewinding: @@ -542,7 +542,8 @@ def solve_equilibrium_problem(self, crve, material_state, mac_load, # Clusters strain initial iterative guess # ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ # Set clusters strain initial iterative guess - global_strain_mf = self._init_global_strain_mf() + global_strain_mf = \ + self._init_global_strain_mf(crve, material_state) # ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ # Set far-field strain initial iterative guess farfield_strain_mf = self._init_farfield_strain_mf() @@ -860,9 +861,9 @@ def solve_equilibrium_problem(self, crve, material_state, mac_load, is_trigger, target_clusters, target_clusters_data = \ adaptivity_manager.get_target_clusters( crve.get_phase_clusters(), crve.get_voxels_clusters(), + material_state.get_clusters_state(), material_state.get_clusters_def_gradient_mf(), material_state.get_clusters_def_gradient_old_mf(), - material_state.get_clusters_state(), material_state.get_clusters_state_old(), clusters_sct_mf, clusters_sct_old_mf, clusters_residuals_mf, inc, @@ -1059,7 +1060,8 @@ def solve_equilibrium_problem(self, crve, material_state, mac_load, # Update clustering adaptivity output file voxels_output.write_voxels_output_file( self._n_dim, comp_order, crve, - material_state.get_clusters_state()) + material_state.get_clusters_state(), + material_state.get_clusters_def_gradient_mf()) # Increment post-processing time self._post_process_time += time.time() - procedure_init_time # @@ -1130,11 +1132,16 @@ def solve_equilibrium_problem(self, crve, material_state, mac_load, # Set increment initial time inc_init_time = time.time() # ------------------------------------------------------------------------- - def _init_global_strain_mf(self, mode='last_converged'): + def _init_global_strain_mf(self, crve, material_state, + mode='last_converged'): """Set clusters strains initial iterative guess. Parameters ---------- + crve : CRVE + Cluster-Reduced Representative Volume Element. + material_state : MaterialState + CRVE material constitutive state at rewind state. mode : {'last_converged',}, default='last_converged' Strategy to set clusters incremental strains initial iterative guess. @@ -1144,10 +1151,58 @@ def _init_global_strain_mf(self, mode='last_converged'): global_strain_mf : numpy.ndarray (1d) Global vector of clusters strains stored in matricial form. """ + # Set strain/stress components order according to problem strain + # formulation + if self._strain_formulation == 'infinitesimal': + comp_order = self._comp_order_sym + elif self._strain_formulation == 'finite': + comp_order = self._comp_order_nsym + else: + raise RuntimeError('Unknown problem strain formulation.') + # ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + # Get material phases + material_phases = material_state.get_material_phases() + # Get clusters associated with each material phase + phase_clusters = crve.get_phase_clusters() + # Get total number of clusters + n_total_clusters = crve.get_n_total_clusters() + # ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + # Set clusters strain initial iterative guess associated with the last + # converged solution if mode == 'last_converged': - # Set initial iterative guess associated with the last converged - # solution - global_strain_mf = copy.deepcopy(self._global_strain_old_mf) + # Initialize initial iterative guess + global_strain_mf = np.zeros((n_total_clusters*len(comp_order))) + # ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + # Get last converged clusters state variables + clusters_state_old = material_state.get_clusters_state_old() + # Get last converged clusters deformation gradient + clusters_def_gradient_old_mf = \ + material_state.get_clusters_def_gradient_old_mf() + # ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + # Initialize material cluster strain range indexes + i_init = 0 + i_end = i_init + len(comp_order) + # Loop over material phases + for mat_phase in material_phases: + # Loop over material phase clusters + for cluster in phase_clusters[mat_phase]: + # Get last converged material cluster infinitesimal strain + # tensor (infinitesimal strains) or deformation gradient + # (finite strains) + if self._strain_formulation == 'infinitesimal': + strain_old_mf = \ + clusters_state_old[str(cluster)]['strain_mf'] + else: + strain_old_mf = \ + clusters_def_gradient_old_mf[str(cluster)] + # ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + # Assemble to initial iterative guess + global_strain_mf[i_init:i_end] = strain_old_mf + # ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + # Update cluster strain range indexes + i_init = i_init + len(comp_order) + i_end = i_init + len(comp_order) + # ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ else: raise RuntimeError('Unavailable strategy to set clusters strains ' 'initial iterative guess.') diff --git a/src/cratepy/tensor/matrixoperations.py b/src/cratepy/tensor/matrixoperations.py index 553d0fa2..2a5d1d75 100644 --- a/src/cratepy/tensor/matrixoperations.py +++ b/src/cratepy/tensor/matrixoperations.py @@ -405,7 +405,7 @@ def kelvin_factor(idx, comp_order): # # Matrix condensation # ============================================================================= -def get_condensed_matrix(matrix, rows, cols): +def get_condensed_matrix(matrix, rows, cols, is_check_indexes=False): """Perform condensation of matrix given a set of rows and columns. Parameters @@ -416,6 +416,10 @@ def get_condensed_matrix(matrix, rows, cols): Indexes of rows to keep in condensed matrix. cols : numpy.ndarray (1d) Indexes of columns to keep in condensed matrix. + is_check_indexes : bool, default=False + If True, then check validity of condensed rows and columns indexes + before performing matrix condensation. May yield non-negligible + overhead cost when condensing large matrix. Returns ------- @@ -423,24 +427,26 @@ def get_condensed_matrix(matrix, rows, cols): Condensed matrix. """ # Check validity of rows and columns indexes to perform the condensation - if not np.all([isinstance(rows[i], int) or isinstance(rows[i], np.integer) - for i in range(len(rows))]): - raise RuntimeError('All the indexes specified to perform a matrix ' - 'condensation must be non-negative integers.') - elif not np.all([isinstance(cols[i], int) - or isinstance(cols[i], np.integer) - for i in range(len(cols))]): - raise RuntimeError('All the indexes specified to perform a matrix ' - 'condensation must be non-negative integers.') - elif len(list(dict.fromkeys(rows))) != len(rows) or \ - len(list(dict.fromkeys(cols))) != len(cols): - raise RuntimeError('Duplicated rows or columns indexes.') - elif np.any([rows[i] not in range(matrix.shape[0]) - for i in range(len(rows))]): - raise RuntimeError('Out-of-bounds row index.') - elif np.any([cols[i] not in range(matrix.shape[1]) - for i in range(len(cols))]): - raise RuntimeError('Out-of-bounds column index.') + if is_check_indexes: + if not np.all([isinstance(rows[i], int) + or isinstance(rows[i], np.integer) + for i in range(len(rows))]): + raise RuntimeError('All the indexes specified to perform a matrix ' + 'condensation must be non-negative integers.') + elif not np.all([isinstance(cols[i], int) + or isinstance(cols[i], np.integer) + for i in range(len(cols))]): + raise RuntimeError('All the indexes specified to perform a matrix ' + 'condensation must be non-negative integers.') + elif len(list(dict.fromkeys(rows))) != len(rows) or \ + len(list(dict.fromkeys(cols))) != len(cols): + raise RuntimeError('Duplicated rows or columns indexes.') + elif np.any([rows[i] not in range(matrix.shape[0]) + for i in range(len(rows))]): + raise RuntimeError('Out-of-bounds row index.') + elif np.any([cols[i] not in range(matrix.shape[1]) + for i in range(len(cols))]): + raise RuntimeError('Out-of-bounds column index.') # ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ # Build auxiliary matrices with rows and columns condensation indexes rows_matrix = np.zeros((len(rows), len(cols)), dtype=int)