diff --git a/.gitattributes b/.gitattributes index 1437f7412..05b249328 100644 --- a/.gitattributes +++ b/.gitattributes @@ -12,8 +12,12 @@ *.pbz2 binary # ignore files from MSC.Marc in language statistics -installation/mods_MarcMentat/20*/* linguist-vendored +installation/mods_MarcMentat/** linguist-vendored src/Marc/include/* linguist-vendored +installation/mods_MarcMentat/apply_DAMASK_modifications.py linguist-vendored=false # ignore reference files for tests in language statistics -python/tests/reference/* linguist-vendored +python/tests/reference/** linguist-vendored + +# ignore deprecated scripts +processing/legacy/** linguist-vendored diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml index 9df1e9330..3228d8db6 100644 --- a/.gitlab-ci.yml +++ b/.gitlab-ci.yml @@ -233,7 +233,7 @@ source_distribution: stage: deploy script: - cd $(mktemp -d) - - $DAMASKROOT/PRIVATE/releasing/deployMe.sh $CI_COMMIT_SHA + - $DAMASKROOT/PRIVATE/releasing/deploy.sh $DAMASKROOT $CI_COMMIT_SHA except: - master - release diff --git a/PRIVATE b/PRIVATE index 1490f9741..afffa8d04 160000 --- a/PRIVATE +++ b/PRIVATE @@ -1 +1 @@ -Subproject commit 1490f97417664d6ae11d7ceafb6b98c9cc2dade1 +Subproject commit afffa8d04e110282e514a4e57d0bad9c76effe01 diff --git a/VERSION b/VERSION index 01ac5403c..543799db9 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -v3.0.0-alpha2-866-g1be1a72a0 +v3.0.0-alpha3-2-gd0dd1fd83 diff --git a/examples/config/phase/mechanical/plastic/dislotwin_IF-steel.yaml b/examples/config/phase/mechanical/plastic/dislotwin_IF-steel.yaml index 9d0c5be8b..bebe0f341 100644 --- a/examples/config/phase/mechanical/plastic/dislotwin_IF-steel.yaml +++ b/examples/config/phase/mechanical/plastic/dislotwin_IF-steel.yaml @@ -5,6 +5,7 @@ references: 10.1016/j.ijplas.2020.102779 - K. Sedighiani et al., Mechanics of Materials, submitted +output: [rho_dip, rho_mob] N_sl: [12, 12] b_sl: [2.49e-10, 2.49e-10] rho_mob_0: [2.81e12, 2.8e12] diff --git a/python/damask/__init__.py b/python/damask/__init__.py index af9933954..ad46d454f 100644 --- a/python/damask/__init__.py +++ b/python/damask/__init__.py @@ -1,11 +1,4 @@ -""" -Tools for pre and post processing of DAMASK simulations. - -Modules that contain only one class (of the same name), -are prefixed by a '_'. For example, '_colormap' contains -a class called 'Colormap' which is imported as 'damask.Colormap'. - -""" +"""Tools for managing DAMASK simulations.""" from pathlib import Path as _Path import re as _re @@ -15,7 +8,6 @@ with open(_Path(__file__).parent/_Path('VERSION')) as _f: version = _re.sub(r'^v','',_f.readline().strip()) __version__ = version -# make classes directly accessible as damask.Class from . import util # noqa from . import seeds # noqa from . import tensor # noqa @@ -23,6 +15,8 @@ from . import mechanics # noqa from . import solver # noqa from . import grid_filters # noqa from . import lattice # noqa +#Modules that contain only one class (of the same name), are prefixed by a '_'. +#For example, '_colormap' containsa class called 'Colormap' which is imported as 'damask.Colormap'. from ._rotation import Rotation # noqa from ._orientation import Orientation # noqa from ._table import Table # noqa diff --git a/python/damask/_colormap.py b/python/damask/_colormap.py index de63508b3..02e2b922c 100644 --- a/python/damask/_colormap.py +++ b/python/damask/_colormap.py @@ -91,6 +91,16 @@ class Colormap(mpl.colors.ListedColormap): - 'lab': CIE Lab. - 'msh': Msh (for perceptual uniform interpolation). + Returns + ------- + new : damask.Colormap + Colormap within given bounds. + + Examples + -------- + >>> import damask + >>> damask.Colormap.from_range((0,0,1),(0,0,0),'blue_to_black') + """ low_high = np.vstack((low,high)) if model.lower() == 'rgb': @@ -150,6 +160,16 @@ class Colormap(mpl.colors.ListedColormap): This parameter is not used for matplotlib colormaps that are of type `ListedColormap`. + Returns + ------- + new : damask.Colormap + Predefined colormap. + + Examples + -------- + >>> import damask + >>> damask.Colormap.from_predefined('strain') + """ # matplotlib presets try: @@ -220,6 +240,11 @@ class Colormap(mpl.colors.ListedColormap): damask.Colormap The reversed colormap. + Examples + -------- + >>> import damask + >>> damask.Colormap.from_predefined('stress').reversed() + """ rev = super(Colormap,self).reversed(name) return Colormap(np.array(rev.colors),rev.name[:-4] if rev.name.endswith('_r_r') else rev.name) @@ -239,8 +264,8 @@ class Colormap(mpl.colors.ListedColormap): Returns ------- - f - File handle + f : file object + File handle with write access. """ if fname is None: diff --git a/python/damask/_configmaterial.py b/python/damask/_configmaterial.py index 557907594..e927874f3 100644 --- a/python/damask/_configmaterial.py +++ b/python/damask/_configmaterial.py @@ -14,7 +14,7 @@ class ConfigMaterial(Config): A complete material configuration file has the entries 'material', 'phase', and 'homogenization'. For use in DAMASK, it needs to be stored as 'material.yaml'. - + """ def __init__(self,d=None): @@ -60,54 +60,6 @@ class ConfigMaterial(Config): return super(ConfigMaterial,cls).load(fname) - @staticmethod - def from_table(table,**kwargs): - """ - Generate from an ASCII table. - - Parameters - ---------- - table : damask.Table - Table that contains material information. - **kwargs - Keyword arguments where the key is the name and the value specifies - the label of the data column in the table. - - Examples - -------- - >>> import damask - >>> import damask.ConfigMaterial as cm - >>> t = damask.Table.load('small.txt') - >>> t - pos pos pos qu qu qu qu phase homog - 0 0 0 0 0.19 0.8 0.24 -0.51 Aluminum SX - 1 1 0 0 0.8 0.19 0.24 -0.51 Steel SX - 1 1 1 0 0.8 0.19 0.24 -0.51 Steel SX - >>> cm.from_table(t,O='qu',phase='phase',homogenization='homog') - material: - - constituents: - - O: [0.19, 0.8, 0.24, -0.51] - v: 1.0 - phase: Aluminum - homogenization: SX - - constituents: - - O: [0.8, 0.19, 0.24, -0.51] - v: 1.0 - phase: Steel - homogenization: SX - homogenization: {} - phase: {} - - """ - kwargs_ = {k:table.get(v) for k,v in kwargs.items()} - - _,idx = np.unique(np.hstack(list(kwargs_.values())),return_index=True,axis=0) - idx = np.sort(idx) - kwargs_ = {k:np.atleast_1d(v[idx].squeeze()) for k,v in kwargs_.items()} - - return ConfigMaterial().material_add(**kwargs_) - - @staticmethod def load_DREAM3D(fname, grain_data=None,cell_data=None,cell_ensemble_data='CellEnsembleData', @@ -181,9 +133,69 @@ class ConfigMaterial(Config): return base_config.material_add(**constituent,homogenization='direct') + @staticmethod + def from_table(table,**kwargs): + """ + Generate from an ASCII table. + + Parameters + ---------- + table : damask.Table + Table that contains material information. + **kwargs + Keyword arguments where the key is the name and the value specifies + the label of the data column in the table. + + Examples + -------- + >>> import damask + >>> import damask.ConfigMaterial as cm + >>> t = damask.Table.load('small.txt') + >>> t + pos pos pos qu qu qu qu phase homog + 0 0 0 0 0.19 0.8 0.24 -0.51 Aluminum SX + 1 1 0 0 0.8 0.19 0.24 -0.51 Steel SX + 1 1 1 0 0.8 0.19 0.24 -0.51 Steel SX + >>> cm.from_table(t,O='qu',phase='phase',homogenization='homog') + material: + - constituents: + - O: [0.19, 0.8, 0.24, -0.51] + v: 1.0 + phase: Aluminum + homogenization: SX + - constituents: + - O: [0.8, 0.19, 0.24, -0.51] + v: 1.0 + phase: Steel + homogenization: SX + homogenization: {} + phase: {} + + """ + kwargs_ = {k:table.get(v) for k,v in kwargs.items()} + + _,idx = np.unique(np.hstack(list(kwargs_.values())),return_index=True,axis=0) + idx = np.sort(idx) + kwargs_ = {k:np.atleast_1d(v[idx].squeeze()) for k,v in kwargs_.items()} + + return ConfigMaterial().material_add(**kwargs_) + + @property def is_complete(self): - """Check for completeness.""" + """ + Check for completeness. + + Only the general file layout is considered. + This check does not consider whether parameters for + a particular phase/homogenization model are missing. + + Returns + ------- + complete : bool + Whether the material.yaml definition is complete. + + """ ok = True for top_level in ['homogenization','phase','material']: ok &= top_level in self @@ -236,7 +248,19 @@ class ConfigMaterial(Config): @property def is_valid(self): - """Check for valid content.""" + """ + Check for valid content. + + Only the generic file content is considered. + This check does not consider whether parameters for a + particular phase/homogenization mode are out of bounds. + + Returns + ------- + valid : bool + Whether the material.yaml definition is valid. + + """ ok = True if 'phase' in self: @@ -282,7 +306,7 @@ class ConfigMaterial(Config): Returns ------- - cfg : damask.ConfigMaterial + updated : damask.ConfigMaterial Updated material configuration. """ @@ -311,7 +335,7 @@ class ConfigMaterial(Config): Returns ------- - cfg : damask.ConfigMaterial + updated : damask.ConfigMaterial Updated material configuration. """ @@ -336,53 +360,61 @@ class ConfigMaterial(Config): Returns ------- - cfg : damask.ConfigMaterial + updated : damask.ConfigMaterial Updated material configuration. Examples -------- + Create a dual-phase steel microstructure for micromechanical simulations: + >>> import numpy as np >>> import damask - >>> m = damask.ConfigMaterial().material_add(phase = ['Aluminum','Steel'], - ... O = damask.Rotation.from_random(2), - ... homogenization = 'SX') + >>> m = damask.ConfigMaterial() + >>> m = m.material_add(phase = ['Ferrite','Martensite'], + ... O = damask.Rotation.from_random(2), + ... homogenization = 'SX') >>> m material: - constituents: - O: [0.577764, -0.146299, -0.617669, 0.513010] v: 1.0 - phase: Aluminum + phase: Ferrite homogenization: SX - constituents: - O: [0.184176, 0.340305, 0.737247, 0.553840] v: 1.0 - phase: Steel + phase: Martensite homogenization: SX homogenization: {} phase: {} - >>> m = damask.ConfigMaterial().material_add(phase = np.array(['Austenite','Martensite']).reshape(1,2), - ... O = damask.Rotation.from_random((2,2)), - ... v = np.array([0.2,0.8]).reshape(1,2), - ... homogenization = ['A','B']) + Create a duplex stainless steel microstructure for forming simulations: + + >>> import numpy as np + >>> import damask + >>> m = damask.ConfigMaterial() + >>> m = m.material_add(phase = np.array(['Austenite','Ferrite']).reshape(1,2), + ... O = damask.Rotation.from_random((2,2)), + ... v = np.array([0.2,0.8]).reshape(1,2), + ... homogenization = 'Taylor') >>> m material: - constituents: - phase: Austenite O: [0.659802978293224, 0.6953785848195171, 0.22426295326327111, -0.17554139512785227] v: 0.2 - - phase: Martensite + - phase: Ferrite O: [0.49356745891301596, 0.2841806579193434, -0.7487679215072818, -0.339085707289975] v: 0.8 - homogenization: A + homogenization: Taylor - constituents: - phase: Austenite O: [0.26542221365204055, 0.7268854930702071, 0.4474726435701472, -0.44828201137283735] v: 0.2 - - phase: Martensite + - phase: Ferrite O: [0.6545817158479885, -0.08004812803625233, -0.6226561293931374, 0.4212059104577611] v: 0.8 - homogenization: B + homogenization: Taylor homogenization: {} phase: {} diff --git a/python/damask/_grid.py b/python/damask/_grid.py index 308be5c80..aebd8f92c 100644 --- a/python/damask/_grid.py +++ b/python/damask/_grid.py @@ -161,6 +161,11 @@ class Grid: Grid file to read. Valid extension is .vtr, which will be appended if not given. + Returns + ------- + loaded : damask.Grid + Grid-based geometry from file. + """ v = VTK.load(fname if str(fname).endswith('.vtr') else str(fname)+'.vtr') comments = v.get_comments() @@ -190,6 +195,11 @@ class Grid: fname : str, pathlib.Path, or file handle Geometry file to read. + Returns + ------- + loaded : damask.Grid + Grid-based geometry from file. + """ warnings.warn('Support for ASCII-based geom format will be removed in DAMASK 3.1.0', DeprecationWarning,2) try: @@ -281,6 +291,10 @@ class Grid: and grain- or cell-wise data. Defaults to None, in which case it is set as the path that contains _SIMPL_GEOMETRY/SPACING. + Returns + ------- + loaded : damask.Grid + Grid-based geometry from file. """ b = util.DREAM3D_base_group(fname) if base_group is None else base_group @@ -319,6 +333,11 @@ class Grid: Label(s) of the columns containing the material definition. Each unique combination of values results in one material ID. + Returns + ------- + new : damask.Grid + Grid-based geometry from values in table. + """ cells,size,origin = grid_filters.cellsSizeOrigin_coordinates0_point(table.get(coordinates)) @@ -356,6 +375,11 @@ class Grid: periodic : Boolean, optional Assume grid to be periodic. Defaults to True. + Returns + ------- + new : damask.Grid + Grid-based geometry from tessellation. + """ if periodic: weights_p = np.tile(weights,27) # Laguerre weights (1,2,3,1,2,3,...,1,2,3) @@ -405,6 +429,11 @@ class Grid: periodic : Boolean, optional Assume grid to be periodic. Defaults to True. + Returns + ------- + new : damask.Grid + Grid-based geometry from tessellation. + """ coords = grid_filters.coordinates0_point(cells,size).reshape(-1,3) KDTree = spatial.cKDTree(seeds,boxsize=size) if periodic else spatial.cKDTree(seeds) @@ -478,6 +507,11 @@ class Grid: materials : (int, int), optional Material IDs. Defaults to (1,2). + Returns + ------- + new : damask.Grid + Grid-based geometry from definition of minimal surface. + Notes ----- The following triply-periodic minimal surfaces are implemented: @@ -598,6 +632,11 @@ class Grid: periodic : Boolean, optional Assume grid to be periodic. Defaults to True. + Returns + ------- + updated : damask.Grid + Updated grid-based geometry. + """ # radius and center r = np.array(dimension)/2.0*self.size/self.cells if np.array(dimension).dtype in np.sctypes['int'] else \ @@ -638,6 +677,11 @@ class Grid: reflect : bool, optional Reflect (include) outermost layers. Defaults to False. + Returns + ------- + updated : damask.Grid + Updated grid-based geometry. + """ valid = ['x','y','z'] if not set(directions).issubset(valid): @@ -670,6 +714,11 @@ class Grid: Direction(s) along which the grid is flipped. Valid entries are 'x', 'y', 'z'. + Returns + ------- + updated : damask.Grid + Updated grid-based geometry. + """ valid = ['x','y','z'] if not set(directions).issubset(valid): @@ -695,6 +744,11 @@ class Grid: periodic : Boolean, optional Assume grid to be periodic. Defaults to True. + Returns + ------- + updated : damask.Grid + Updated grid-based geometry. + """ return Grid(material = ndimage.interpolation.zoom( self.material, @@ -723,6 +777,11 @@ class Grid: periodic : Boolean, optional Assume grid to be periodic. Defaults to True. + Returns + ------- + updated : damask.Grid + Updated grid-based geometry. + """ def mostFrequent(arr,selection=None): me = arr[arr.size//2] @@ -746,7 +805,15 @@ class Grid: def renumber(self): - """Renumber sorted material indices as 0,...,N-1.""" + """ + Renumber sorted material indices as 0,...,N-1. + + Returns + ------- + updated : damask.Grid + Updated grid-based geometry. + + """ _,renumbered = np.unique(self.material,return_inverse=True) return Grid(material = renumbered.reshape(self.cells), @@ -767,6 +834,11 @@ class Grid: fill : int or float, optional Material index to fill the corners. Defaults to material.max() + 1. + Returns + ------- + updated : damask.Grid + Updated grid-based geometry. + """ if fill is None: fill = np.nanmax(self.material) + 1 dtype = float if isinstance(fill,float) or self.material.dtype in np.sctypes['float'] else int @@ -802,6 +874,11 @@ class Grid: fill : int or float, optional Material index to fill the background. Defaults to material.max() + 1. + Returns + ------- + updated : damask.Grid + Updated grid-based geometry. + """ if offset is None: offset = 0 if fill is None: fill = np.nanmax(self.material) + 1 @@ -834,6 +911,11 @@ class Grid: to_material : iterable of ints New material indices. + Returns + ------- + updated : damask.Grid + Updated grid-based geometry. + """ def mp(entry,mapper): return mapper[entry] if entry in mapper else entry @@ -849,7 +931,15 @@ class Grid: def sort(self): - """Sort material indices such that min(material) is located at (0,0,0).""" + """ + Sort material indices such that min(material) is located at (0,0,0). + + Returns + ------- + updated : damask.Grid + Updated grid-based geometry. + + """ a = self.material.flatten(order='F') from_ma = pd.unique(a) sort_idx = np.argsort(from_ma) @@ -884,6 +974,11 @@ class Grid: periodic : Boolean, optional Assume grid to be periodic. Defaults to True. + Returns + ------- + updated : damask.Grid + Updated grid-based geometry. + """ def tainted_neighborhood(stencil,trigger): diff --git a/python/damask/_orientation.py b/python/damask/_orientation.py index 5cc412c98..00e74d36c 100644 --- a/python/damask/_orientation.py +++ b/python/damask/_orientation.py @@ -8,15 +8,15 @@ from . import tensor _parameter_doc = \ """lattice : str - Either a crystal family out of [triclinic, monoclinic, orthorhombic, tetragonal, hexagonal, cubic] - or a Bravais lattice out of [aP, mP, mS, oP, oS, oI, oF, tP, tI, hP, cP, cI, cF]. - When specifying a Bravais lattice, additional lattice parameters might be required: + Either a crystal family out of {triclinic, monoclinic, orthorhombic, tetragonal, hexagonal, cubic} + or a Bravais lattice out of {aP, mP, mS, oP, oS, oI, oF, tP, tI, hP, cP, cI, cF}. + When specifying a Bravais lattice, additional lattice parameters might be required. a : float, optional - Length of lattice parameter "a". + Length of lattice parameter 'a'. b : float, optional - Length of lattice parameter "b". + Length of lattice parameter 'b'. c : float, optional - Length of lattice parameter "c". + Length of lattice parameter 'c'. alpha : float, optional Angle between b and c lattice basis. beta : float, optional diff --git a/python/damask/_result.py b/python/damask/_result.py index 86c488416..2c6cb7285 100644 --- a/python/damask/_result.py +++ b/python/damask/_result.py @@ -57,11 +57,30 @@ def _empty_like(dataset,N_materialpoints,fill_float,fill_int): class Result: """ - Manipulate and read DADF5 files. + Add data to and export data from a DADF5 file. + + A DADF5 (DAMASK HDF5) file contain DAMASK results. + Its group/folder structure reflects the layout in material.yaml. + + This class provides a customable view on the DADF5 file. + Upon initialization, all attributes are visible. + Derived quantities are added to the file and existing data is + exported based on the current view. + + Examples + -------- + Open 'my_file.hdf5', which needs to contain deformation gradient 'F' + and first Piola-Kirchhoff stress 'P', add the Mises equivalent of the + Cauchy stress, and export it to VTK (file) and numpy.ndarray (memory). + + >>> import damask + >>> r = damask.Result('my_file.hdf5') + >>> r.add_Cauchy() + >>> r.add_equivalent_Mises('sigma') + >>> r.save_VTK() + >>> r_last = r.view('increments',-1) + >>> sigma_vM_last = r_last.get('sigma_vM') - DADF5 (DAMASK HDF5) files contain DAMASK results. - The group/folder structure reflects the input data - in material.yaml. """ def __init__(self,fname): @@ -95,10 +114,10 @@ class Result: self.N_materialpoints, self.N_constituents = np.shape(f[f'cell_to/phase']) - self.homogenizations = [m.decode() for m in np.unique(f[f'cell_to/homogenization']['label'])] - self.homogenizations = sorted(self.homogenizations,key=util.natural_sort) - self.phases = [c.decode() for c in np.unique(f[f'cell_to/phase']['label'])] - self.phases = sorted(self.phases,key=util.natural_sort) + self.homogenization = f[f'cell_to/homogenization']['label'].astype('str') + self.homogenizations = sorted(np.unique(self.homogenization),key=util.natural_sort) + self.phase = f[f'cell_to/phase']['label'].astype('str') + self.phases = sorted(np.unique(self.phase),key=util.natural_sort) self.fields = [] for c in self.phases: @@ -154,6 +173,11 @@ class Result: Name of datasets; supports '?' and '*' wildcards. True is equivalent to '*', False is equivalent to []. + Returns + ------- + view : damask.Result + Modified or new view on the DADF5 file. + """ # allow True/False and string arguments if datasets is True: @@ -166,6 +190,7 @@ class Result: if what == 'increments': choice = [c if isinstance(c,str) and c.startswith('increment_') else f'increment_{c}' for c in choice] + if datasets == -1: choice = [self.increments[-1]] elif what == 'times': what = 'increments' if choice == ['*']: @@ -197,15 +222,31 @@ class Result: return dup - def allow_modification(self): - """Allow to overwrite existing data.""" + def modification_enable(self): + """ + Allow to modify existing data. + + Returns + ------- + modified_view : damask.Result + View where data is not write-protected. + + """ print(util.warn('Warning: Modification of existing datasets allowed!')) dup = self.copy() dup._allow_modification = True return dup - def disallow_modification(self): - """Disallow to overwrite existing data (default case).""" + def modification_disable(self): + """ + Disallow to modify existing data (default case). + + Returns + ------- + modified_view : damask.Result + View where data is write-protected. + + """ dup = self.copy() dup._allow_modification = False return dup @@ -213,7 +254,7 @@ class Result: def increments_in_range(self,start,end): """ - Select all increments within a given range. + Get all increments within a given range. Parameters ---------- @@ -222,6 +263,11 @@ class Result: end : int or str End increment. + Returns + ------- + increments : list of ints + Increment number of all increments within the given bounds. + """ selected = [] for i,inc in enumerate([int(i[10:]) for i in self.increments]): @@ -233,7 +279,7 @@ class Result: def times_in_range(self,start,end): """ - Select all increments within a given time range. + Get all increments within a given time range. Parameters ---------- @@ -242,6 +288,11 @@ class Result: end : float Time of end increment. + Returns + ------- + times : list of float + Simulation time of all increments within the given bounds. + """ selected = [] for i,time in enumerate(self.times): @@ -262,6 +313,25 @@ class Result: Name of datasets; supports '?' and '*' wildcards. True is equivalent to '*', False is equivalent to []. + Returns + ------- + view : damask.Result + View with where selected attributes are visible. + + Examples + -------- + Get a view that shows only results from the initial configuration: + + >>> import damask + >>> r = damask.Result('my_file.hdf5') + >>> r_first = r.view('increment',0) + + Get a view that shows all results of in simulation time [10,40]: + + >>> import damask + >>> r = damask.Result('my_file.hdf5') + >>> r_t10to40 = r.view('times',r.times_in_range(10.0,40.0)) + """ return self._manage_view('set',what,datasets) @@ -278,6 +348,20 @@ class Result: Name of datasets; supports '?' and '*' wildcards. True is equivalent to '*', False is equivalent to []. + Returns + ------- + modified_view : damask.Result + View with more visible attributes. + + Examples + -------- + Get a view that shows only results from first and last increment: + + >>> import damask + >>> r_empty = damask.Result('my_file.hdf5').view('increments',False) + >>> r_first = r_empty.view_more('increments',0) + >>> r_first_and_last = r.first.view_more('increments',-1) + """ return self._manage_view('add',what,datasets) @@ -294,37 +378,98 @@ class Result: Name of datasets; supports '?' and '*' wildcards. True is equivalent to '*', False is equivalent to []. + Returns + ------- + modified_view : damask.Result + View with less visible attributes. + + Examples + -------- + Get a view that does not show the undeformed configuration: + + >>> import damask + >>> r_all = damask.Result('my_file.hdf5') + >>> r_deformed = r_all.view_less('increments',0) + """ return self._manage_view('del',what,datasets) - def rename(self,name_old,name_new): + def rename(self,name_src,name_dst): """ - Rename dataset. + Rename/move datasets (within the same group/folder). + + This operation is discouraged because the history of the + data becomes untracable and scientific integrity cannot be + ensured. Parameters ---------- - name_old : str - Name of the dataset to be renamed. - name_new : str - New name of the dataset. + name_src : str + Name of the datasets to be renamed. + name_dst : str + New name of the datasets. + + Examples + -------- + Rename datasets containing the deformation gradient from 'F' to 'def_grad': + + >>> import damask + >>> r = damask.Result('my_file.hdf5') + >>> r_unprotected = r.modification_enable() + >>> r_unprotected.rename('F','def_grad') """ if not self._allow_modification: - raise PermissionError('Rename operation not permitted') + raise PermissionError('Renaming datasets not permitted') with h5py.File(self.fname,'a') as f: for inc in self.visible['increments']: for ty in ['phase','homogenization']: for label in self.visible[ty+'s']: - for field in self.visible['fields']: - path_old = '/'.join([inc,ty,label,field,name_old]) - path_new = '/'.join([inc,ty,label,field,name_new]) - if path_old in f.keys(): - f[path_new] = f[path_old] - f[path_new].attrs['renamed'] = f'original name: {name_old}' if h5py3 else \ - f'original name: {name_old}'.encode() - del f[path_old] + for field in _match(self.visible['fields'],f['/'.join([inc,ty,label])].keys()): + path_src = '/'.join([inc,ty,label,field,name_src]) + path_dst = '/'.join([inc,ty,label,field,name_dst]) + if path_src in f.keys(): + f[path_dst] = f[path_src] + f[path_dst].attrs['renamed'] = f'original name: {name_src}' if h5py3 else \ + f'original name: {name_src}'.encode() + del f[path_src] + + + def remove(self,name): + """ + Remove/delete datasets. + + This operation is discouraged because the history of the + data becomes untracable and scientific integrity cannot be + ensured. + + Parameters + ---------- + name : str + Name of the datasets to be deleted. + + Examples + -------- + Delete the deformation gradient 'F': + + >>> import damask + >>> r = damask.Result('my_file.hdf5') + >>> r_unprotected = r.modification_enable() + >>> r_unprotected.remove('F') + + """ + if not self._allow_modification: + raise PermissionError('Removing datasets not permitted') + + with h5py.File(self.fname,'a') as f: + for inc in self.visible['increments']: + for ty in ['phase','homogenization']: + for label in self.visible[ty+'s']: + for field in _match(self.visible['fields'],f['/'.join([inc,ty,label])].keys()): + path = '/'.join([inc,ty,label,field,name]) + if path in f.keys(): del f[path] def list_data(self): @@ -337,7 +482,7 @@ class Result: msg = ' '.join([msg,f'{ty}\n']) for label in self.visible[ty+'s']: msg = ' '.join([msg,f'{label}\n']) - for field in self.visible['fields']: + for field in _match(self.visible['fields'],f['/'.join([inc,ty,label])].keys()): msg = ' '.join([msg,f'{field}\n']) for d in f['/'.join([inc,ty,label,field])].keys(): dataset = f['/'.join([inc,ty,label,field,d])] @@ -357,7 +502,7 @@ class Result: @property def coordinates0_point(self): - """Return initial coordinates of the cell centers.""" + """Initial/undeformed cell center coordinates.""" if self.structured: return grid_filters.coordinates0_point(self.cells,self.size,self.origin).reshape(-1,3,order='F') else: @@ -366,7 +511,7 @@ class Result: @property def coordinates0_node(self): - """Return initial coordinates of the cell centers.""" + """Initial/undeformed nodal coordinates.""" if self.structured: return grid_filters.coordinates0_node(self.cells,self.size,self.origin).reshape(-1,3,order='F') else: @@ -375,6 +520,7 @@ class Result: @property def geometry0(self): + """Initial/undeformed geometry.""" if self.structured: return VTK.from_rectilinear_grid(self.cells,self.size,self.origin) else: @@ -403,7 +549,7 @@ class Result: Parameters ---------- x : str - Label of scalar, vector, or tensor dataset to take absolute value of. + Name of scalar, vector, or tensor dataset to take absolute value of. """ self._add_generic_pointwise(self._add_absolute,{'x':x}) @@ -424,24 +570,51 @@ class Result: 'creator': 'add_calculation' } } - def add_calculation(self,label,formula,unit='n/a',description=None): + def add_calculation(self,name,formula,unit='n/a',description=None): """ Add result of a general formula. Parameters ---------- - label : str - Label of resulting dataset. + name : str + Name of resulting dataset. formula : str - Formula to calculate resulting dataset. Existing datasets are referenced by '#TheirLabel#'. + Formula to calculate resulting dataset. Existing datasets are referenced by '#TheirName#'. unit : str, optional Physical unit of the result. description : str, optional Human-readable description of the result. + Examples + -------- + Add total dislocation density, i.e. the sum of mobile dislocation + density 'rho_mob' and dislocation dipole density 'rho_dip' over + all slip systems: + + >>> import damask + >>> r = damask.Result('my_file.hdf5') + >>> r.add_calculation('rho_mob_total','np.sum(#rho_mob#,axis=1)', + ... '1/m²','total mobile dislocation density') + >>> r.add_calculation('rho_dip_total','np.sum(#rho_dip#,axis=1)', + ... '1/m²','total dislocation dipole density') + >>> r.add_calculation('rho_total','#rho_dip_total#+#rho_mob_total', + ... '1/m²','total dislocation density') + + Add Mises equivalent of the Cauchy stress without storage of + intermediate results. Define a user function for better readability: + + >>> import damask + >>> def equivalent_stress(F,P): + ... sigma = damask.mechanics.stress_Cauchy(F=F,P=P) + ... return damask.mechanics.equivalent_stress_Mises(sigma) + >>> r = damask.Result('my_file.hdf5') + >>> r.enable_user_function(equivalent_stress) + >>> r.add_calculation('sigma_vM','equivalent_stress(#F#,#P#)','Pa', + ... 'Mises equivalent of the Cauchy stress') + """ dataset_mapping = {d:d for d in set(re.findall(r'#(.*?)#',formula))} # datasets used in the formula - args = {'formula':formula,'label':label,'unit':unit,'description':description} + args = {'formula':formula,'label':name,'unit':unit,'description':description} self._add_generic_pointwise(self._add_calculation,dataset_mapping,args) @@ -465,9 +638,9 @@ class Result: Parameters ---------- P : str, optional - Label of the dataset containing the first Piola-Kirchhoff stress. Defaults to 'P'. + Name of the dataset containing the first Piola-Kirchhoff stress. Defaults to 'P'. F : str, optional - Label of the dataset containing the deformation gradient. Defaults to 'F'. + Name of the dataset containing the deformation gradient. Defaults to 'F'. """ self._add_generic_pointwise(self._add_stress_Cauchy,{'P':P,'F':F}) @@ -491,7 +664,15 @@ class Result: Parameters ---------- T : str - Label of tensor dataset. + Name of tensor dataset. + + Examples + -------- + Add the determinant of plastic deformation gradient 'F_p': + + >>> import damask + >>> r = damask.Result('my_file.hdf5') + >>> r.add_determinant('F_p') """ self._add_generic_pointwise(self._add_determinant,{'T':T}) @@ -515,7 +696,7 @@ class Result: Parameters ---------- T : str - Label of tensor dataset. + Name of tensor dataset. """ self._add_generic_pointwise(self._add_deviator,{'T':T}) @@ -546,7 +727,7 @@ class Result: Parameters ---------- T_sym : str - Label of symmetric tensor dataset. + Name of symmetric tensor dataset. eigenvalue : str, optional Eigenvalue. Select from 'max', 'mid', 'min'. Defaults to 'max'. @@ -579,7 +760,7 @@ class Result: Parameters ---------- T_sym : str - Label of symmetric tensor dataset. + Name of symmetric tensor dataset. eigenvalue : str, optional Eigenvalue to which the eigenvector corresponds. Select from 'max', 'mid', 'min'. Defaults to 'max'. @@ -619,9 +800,17 @@ class Result: l : numpy.array of shape (3) Lab frame direction for inverse pole figure. q : str - Label of the dataset containing the crystallographic orientation as quaternions. + Name of the dataset containing the crystallographic orientation as quaternions. Defaults to 'O'. + Examples + -------- + Add the IPF color along [0,1,1] for orientation 'O': + + >>> import damask + >>> r = damask.Result('my_file.hdf5') + >>> r.add_IPF_color(np.array([0,1,1])) + """ self._add_generic_pointwise(self._add_IPF_color,{'q':q},{'l':l}) @@ -644,7 +833,7 @@ class Result: Parameters ---------- T_sym : str - Label of symmetric tensor dataset. + Name of symmetric tensor dataset. """ self._add_generic_pointwise(self._add_maximum_shear,{'T_sym':T_sym}) @@ -678,11 +867,25 @@ class Result: Parameters ---------- T_sym : str - Label of symmetric tensorial stress or strain dataset. + Name of symmetric tensorial stress or strain dataset. kind : {'stress', 'strain', None}, optional Kind of the von Mises equivalent. Defaults to None, in which case it is selected based on the unit of the dataset ('1' -> strain, 'Pa' -> stress). + Examples + -------- + Add the Mises equivalent of the Cauchy stress 'sigma': + + >>> import damask + >>> r = damask.Result('my_file.hdf5') + >>> r.add_equivalent_Mises('sigma') + + Add the Mises equivalent of the spatial logarithmic strain 'epsilon_V^0.0(F)': + + >>> import damask + >>> r = damask.Result('my_file.hdf5') + >>> r.add_equivalent_Mises('epsilon_V^0.0(F)') + """ self._add_generic_pointwise(self._add_equivalent_Mises,{'T_sym':T_sym},{'kind':kind}) @@ -717,7 +920,7 @@ class Result: Parameters ---------- x : str - Label of vector or tensor dataset. + Name of vector or tensor dataset. ord : {non-zero int, inf, -inf, 'fro', 'nuc'}, optional Order of the norm. inf means NumPy’s inf object. For details refer to numpy.linalg.norm. @@ -745,9 +948,9 @@ class Result: Parameters ---------- P : str, optional - Label of first Piola-Kirchhoff stress dataset. Defaults to 'P'. + Name of first Piola-Kirchhoff stress dataset. Defaults to 'P'. F : str, optional - Label of deformation gradient dataset. Defaults to 'F'. + Name of deformation gradient dataset. Defaults to 'F'. """ self._add_generic_pointwise(self._add_stress_second_Piola_Kirchhoff,{'P':P,'F':F}) @@ -786,7 +989,7 @@ class Result: # Parameters # ---------- # q : str - # Label of the dataset containing the crystallographic orientation as quaternions. + # Name of the dataset containing the crystallographic orientation as quaternions. # p : numpy.array of shape (3) # Crystallographic direction or plane. # polar : bool, optional @@ -813,8 +1016,16 @@ class Result: Parameters ---------- - F : str, optional - Label of deformation gradient dataset. + F : str + Name of deformation gradient dataset. + + Examples + -------- + Add the rotational part of deformation gradient 'F': + + >>> import damask + >>> r = damask.Result('my_file.hdf5') + >>> r.add_rotation('F') """ self._add_generic_pointwise(self._add_rotation,{'F':F}) @@ -838,7 +1049,15 @@ class Result: Parameters ---------- T : str - Label of tensor dataset. + Name of tensor dataset. + + Examples + -------- + Add the hydrostatic part of the Cauchy stress 'sigma': + + >>> import damask + >>> r = damask.Result('my_file.hdf5') + >>> r.add_spherical('sigma') """ self._add_generic_pointwise(self._add_spherical,{'T':T}) @@ -864,13 +1083,28 @@ class Result: Parameters ---------- F : str, optional - Label of deformation gradient dataset. Defaults to 'F'. + Name of deformation gradient dataset. Defaults to 'F'. t : {'V', 'U'}, optional Type of the polar decomposition, 'V' for left stretch tensor and 'U' for right stretch tensor. Defaults to 'V'. m : float, optional Order of the strain calculation. Defaults to 0.0. + Examples + -------- + Add the Biot strain based on the deformation gradient 'F': + + >>> import damask + >>> r = damask.Result('my_file.hdf5') + >>> r.strain(t='U',m=0.5) + + Add the plastic Euler-Almansi strain based on the + plastic deformation gradient 'F_p': + + >>> import damask + >>> r = damask.Result('my_file.hdf5') + >>> r.strain('F_p','V',-1) + """ self._add_generic_pointwise(self._add_strain,{'F':F},{'t':t,'m':m}) @@ -894,7 +1128,7 @@ class Result: Parameters ---------- F : str, optional - Label of deformation gradient dataset. Defaults to 'F'. + Name of deformation gradient dataset. Defaults to 'F'. t : {'V', 'U'}, optional Type of the polar decomposition, 'V' for left stretch tensor and 'U' for right stretch tensor. Defaults to 'V'. @@ -947,7 +1181,7 @@ class Result: for inc in self.visible['increments']: for ty in ['phase','homogenization']: for label in self.visible[ty+'s']: - for field in self.visible['fields']: + for field in _match(self.visible['fields'],f['/'.join([inc,ty,label])].keys()): group = '/'.join([inc,ty,label,field]) if set(datasets.values()).issubset(f[group].keys()): groups.append(group) @@ -1001,10 +1235,14 @@ class Result: """ Write XDMF file to directly visualize data in DADF5 file. + The XDMF format is only supported for structured grids + with single phase and single constituent. + For other cases use `save_VTK`. + Parameters ---------- output : (list of) str - Labels of the datasets to read. + Names of the datasets included in the XDMF file. Defaults to '*', in which case all datasets are considered. """ @@ -1082,7 +1320,7 @@ class Result: for ty in ['phase','homogenization']: for label in self.visible[ty+'s']: - for field in self.visible['fields']: + for field in _match(self.visible['fields'],f['/'.join([inc,ty,label])].keys()): for out in _match(output,f['/'.join([inc,ty,label,field])].keys()): name = '/'.join([inc,ty,label,field,out]) shape = f[name].shape[1:] @@ -1131,10 +1369,16 @@ class Result: """ Export to VTK cell/point data. + One VTK file per visible increment is created. + For cell data, the VTK format is a rectilinear grid (.vtr) for + grid-based simulations and an unstructured grid (.vtu) for + mesh-baed simulations. For point data, the VTK format is poly + data (.vtp). + Parameters ---------- output : (list of) str, optional - Labels of the datasets to place. + Names of the datasets included in the VTK file. Defaults to '*', in which case all datasets are exported. mode : {'cell', 'point'} Export in cell format or point format. @@ -1212,7 +1456,7 @@ class Result: Parameters ---------- output : (list of) str - Labels of the datasets to read. + Names of the datasets to read. Defaults to '*', in which case all datasets are read. flatten : bool Remove singular levels of the folder hierarchy. @@ -1264,7 +1508,7 @@ class Result: Parameters ---------- output : (list of) str, optional - Labels of the datasets to place. + Names of the datasets to read. Defaults to '*', in which case all datasets are placed. flatten : bool Remove singular levels of the folder hierarchy. diff --git a/python/damask/_rotation.py b/python/damask/_rotation.py index fff53fc73..6b2869033 100644 --- a/python/damask/_rotation.py +++ b/python/damask/_rotation.py @@ -24,7 +24,7 @@ class Rotation: - Euler angle triplets are implemented using the Bunge convention, with angular ranges of [0,2π], [0,π], [0,2π]. - The rotation angle ω is limited to the interval [0,π]. - - The real part of a quaternion is positive, Re(q) > 0 + - The real part of a quaternion is positive, Re(q) ≥ 0 - P = -1 (as default). Examples @@ -357,7 +357,7 @@ class Rotation: Parameters ---------- - other : Rotation or list of Rotations. + other : damask.Rotation """ return self.copy(rotation=np.vstack(tuple(map(lambda x:x.quaternion, @@ -365,12 +365,28 @@ class Rotation: def flatten(self,order = 'C'): - """Flatten array.""" + """ + Flatten array. + + Returns + ------- + flattened : damask.Rotation + Rotation flattened to single dimension. + + """ return self.copy(rotation=self.quaternion.reshape((-1,4),order=order)) def reshape(self,shape,order = 'C'): - """Reshape array.""" + """ + Reshape array. + + Returns + ------- + reshaped : damask.Rotation + Rotation of given shape. + + """ if isinstance(shape,(int,np.integer)): shape = (shape,) return self.copy(rotation=self.quaternion.reshape(tuple(shape)+(4,),order=order)) @@ -387,6 +403,11 @@ class Rotation: Where to preferentially locate missing dimensions. Either 'left' or 'right' (default). + Returns + ------- + broadcasted : damask.Rotation + Rotation broadcasted to given shape. + """ if isinstance(shape,(int,np.integer)): shape = (shape,) return self.copy(rotation=np.broadcast_to(self.quaternion.reshape(util.shapeshifter(self.shape,shape,mode)+(4,)), @@ -404,7 +425,7 @@ class Rotation: Returns ------- - average : Rotation + average : damask.Rotation Weighted average of original Rotation field. References @@ -438,9 +459,14 @@ class Rotation: Parameters ---------- - other : Rotation + other : damask.Rotation Rotation to which the misorientation is computed. + Returns + ------- + g : damask.Rotation + Misorientation. + """ return other*~self @@ -531,10 +557,10 @@ class Rotation: Returns ------- - rho : numpy.ndarray of shape (...,4) containing + rho : numpy.ndarray of shape (...,4) containing [n_1, n_2, n_3, tan(ω/2)], ǀnǀ = 1 and ω ∈ [0,π] unless compact == True: - numpy.ndarray of shape (...,3) containing + numpy.ndarray of shape (...,3) containing tan(ω/2) [n_1, n_2, n_3], ω ∈ [0,π]. """ diff --git a/python/damask/_table.py b/python/damask/_table.py index 01d8c9863..dd6a981b7 100644 --- a/python/damask/_table.py +++ b/python/damask/_table.py @@ -49,7 +49,7 @@ class Table: Returns ------- - slice : Table + slice : damask.Table Sliced part of the Table. Examples @@ -157,7 +157,7 @@ class Table: Parameters ---------- - other : Table + other : damask.Table Table to compare against. rtol : float, optional Relative tolerance of equality. @@ -185,7 +185,7 @@ class Table: Parameters ---------- - other : Table + other : damask.Table Table to compare against. rtol : float, optional Relative tolerance of equality. @@ -223,6 +223,11 @@ class Table: fname : file, str, or pathlib.Path Filename or file for reading. + Returns + ------- + loaded : damask.Table + Table data from file. + """ try: f = open(fname) @@ -275,6 +280,11 @@ class Table: fname : file, str, or pathlib.Path Filename or file for reading. + Returns + ------- + loaded : damask.Table + Table data from file. + """ try: f = open(fname) @@ -334,14 +344,14 @@ class Table: ---------- label : str Column label. - data : np.ndarray + data : numpy.ndarray New data. info : str, optional Human-readable information about the new data. Returns ------- - table : Table + updated : damask.Table Updated table. """ @@ -367,14 +377,14 @@ class Table: ---------- label : str Column label. - data : np.ndarray + data : numpy.ndarray Modified data. info : str, optional Human-readable information about the modified data. Returns ------- - table : Table + udated : damask.Table Updated table. """ @@ -402,7 +412,7 @@ class Table: Returns ------- - table : Table + udated : damask.Table Updated table. """ @@ -425,7 +435,7 @@ class Table: Returns ------- - table : Table + udated : damask.Table Updated table. """ @@ -451,7 +461,7 @@ class Table: Returns ------- - table : Table + udated : damask.Table Updated table. """ @@ -479,13 +489,13 @@ class Table: Parameters ---------- - other : Table + other : damask.Table Table to append. Returns ------- - table : Table - Concatenated table. + udated : damask.Table + Updated table. """ if self.shapes != other.shapes or not self.data.columns.equals(other.data.columns): @@ -504,13 +514,13 @@ class Table: Parameters ---------- - other : Table + other : damask.Table Table to join. Returns ------- - table : Table - Joined table. + udated : damask.Table + Updated table. """ if set(self.shapes) & set(other.shapes) or self.data.shape[0] != other.data.shape[0]: diff --git a/python/damask/_vtk.py b/python/damask/_vtk.py index 340bd1ea3..2a6cbd3f2 100644 --- a/python/damask/_vtk.py +++ b/python/damask/_vtk.py @@ -3,7 +3,6 @@ import multiprocessing as mp from pathlib import Path import numpy as np -import numpy.ma as ma import vtk from vtk.util.numpy_support import numpy_to_vtk as np_to_vtk from vtk.util.numpy_support import numpy_to_vtkIdTypeArray as np_to_vtkIdTypeArray @@ -51,6 +50,11 @@ class VTK: origin : iterable of float, len (3), optional Spatial origin coordinates. + Returns + ------- + new : damask.VTK + VTK-based geometry without nodal or cell data. + """ vtk_data = vtk.vtkRectilinearGrid() vtk_data.SetDimensions(*(np.array(grid)+1)) @@ -68,7 +72,7 @@ class VTK: """ Create VTK of type vtk.vtkUnstructuredGrid. - This is the common type for FEM solver results. + This is the common type for mesh solver results. Parameters ---------- @@ -80,6 +84,11 @@ class VTK: cell_type : str Name of the vtk.vtkCell subclass. Tested for TRIANGLE, QUAD, TETRA, and HEXAHEDRON. + Returns + ------- + new : damask.VTK + VTK-based geometry without nodal or cell data. + """ vtk_nodes = vtk.vtkPoints() vtk_nodes.SetData(np_to_vtk(nodes)) @@ -110,6 +119,11 @@ class VTK: points : numpy.ndarray of shape (:,3) Spatial position of the points. + Returns + ------- + new : damask.VTK + VTK-based geometry without nodal or cell data. + """ N = points.shape[0] vtk_points = vtk.vtkPoints() @@ -137,8 +151,13 @@ class VTK: fname : str or pathlib.Path Filename for reading. Valid extensions are .vtr, .vtu, .vtp, and .vtk. dataset_type : str, optional - Name of the vtk.vtkDataSet subclass when opening a .vtk file. Valid types are vtkRectilinearGrid, - vtkUnstructuredGrid, and vtkPolyData. + Name of the vtk.vtkDataSet subclass when opening a .vtk file. + Valid types are vtkRectilinearGrid, vtkUnstructuredGrid, and vtkPolyData. + + Returns + ------- + loaded : damask.VTK + VTK-based geometry from file. """ if not os.path.isfile(fname): # vtk has a strange error handling @@ -149,13 +168,13 @@ class VTK: reader.SetFileName(str(fname)) if dataset_type is None: raise TypeError('Dataset type for *.vtk file not given.') - elif dataset_type.lower().endswith('rectilineargrid'): + elif dataset_type.lower().endswith(('rectilineargrid','rectilinear_grid')): reader.Update() vtk_data = reader.GetRectilinearGridOutput() - elif dataset_type.lower().endswith('unstructuredgrid'): + elif dataset_type.lower().endswith(('unstructuredgrid','unstructured_grid')): reader.Update() vtk_data = reader.GetUnstructuredGridOutput() - elif dataset_type.lower().endswith('polydata'): + elif dataset_type.lower().endswith(('polydata','poly_data')): reader.Update() vtk_data = reader.GetPolyDataOutput() else: @@ -246,7 +265,7 @@ class VTK: raise ValueError('No label defined for numpy.ndarray') N_data = data.shape[0] - data_ = np.where(data.mask,data.fill_value,data) if isinstance(data,ma.MaskedArray) else\ + data_ = np.where(data.mask,data.fill_value,data) if isinstance(data,np.ma.MaskedArray) else\ data d = np_to_vtk((data_.astype(np.single) if data_.dtype in [np.double, np.longdouble] else data_).reshape(N_data,-1),deep=True) # avoid large files @@ -277,6 +296,11 @@ class VTK: label : str Data label. + Returns + ------- + data : numpy.ndarray + Data stored under the given label. + """ cell_data = self.vtk_data.GetCellData() for a in range(cell_data.GetNumberOfArrays()): diff --git a/python/damask/grid_filters.py b/python/damask/grid_filters.py index ee929b3bf..bc9194947 100644 --- a/python/damask/grid_filters.py +++ b/python/damask/grid_filters.py @@ -1,15 +1,14 @@ """ Filters for operations on regular grids. -Notes ------ The grids are defined as (x,y,z,...) where x is fastest and z is slowest. This convention is consistent with the layout in grid vtr files. + When converting to/from a plain list (e.g. storage in ASCII table), the following operations are required for tensorial data: -D3 = D1.reshape(cells+(-1,),order='F').reshape(cells+(3,3)) -D1 = D3.reshape(cells+(-1,)).reshape(-1,9,order='F') + - D3 = D1.reshape(cells+(-1,),order='F').reshape(cells+(3,3)) + - D1 = D3.reshape(cells+(-1,)).reshape(-1,9,order='F') """ from scipy import spatial as _spatial @@ -29,10 +28,12 @@ def _ks(size,cells,first_order=False): Correction for first order derivatives, defaults to False. """ - k_sk = _np.where(_np.arange(cells[0])>cells[0]//2,_np.arange(cells[0])-cells[0],_np.arange(cells[0]))/size[0] + k_sk = _np.where(_np.arange(cells[0])>cells[0]//2, + _np.arange(cells[0])-cells[0],_np.arange(cells[0]))/size[0] if cells[0]%2 == 0 and first_order: k_sk[cells[0]//2] = 0 # Nyquist freq=0 for even cells (Johnson, MIT, 2011) - k_sj = _np.where(_np.arange(cells[1])>cells[1]//2,_np.arange(cells[1])-cells[1],_np.arange(cells[1]))/size[1] + k_sj = _np.where(_np.arange(cells[1])>cells[1]//2, + _np.arange(cells[1])-cells[1],_np.arange(cells[1]))/size[1] if cells[1]%2 == 0 and first_order: k_sj[cells[1]//2] = 0 # Nyquist freq=0 for even cells (Johnson, MIT, 2011) k_si = _np.arange(cells[2]//2+1)/size[2] @@ -40,74 +41,89 @@ def _ks(size,cells,first_order=False): return _np.stack(_np.meshgrid(k_sk,k_sj,k_si,indexing = 'ij'), axis=-1) -def curl(size,field): - """ +def curl(size,f): + u""" Calculate curl of a vector or tensor field in Fourier space. Parameters ---------- size : numpy.ndarray of shape (3) Physical size of the periodic field. - field : numpy.ndarray of shape (:,:,:,3) or (:,:,:,3,3) + f : numpy.ndarray of shape (:,:,:,3) or (:,:,:,3,3) Periodic field of which the curl is calculated. + Returns + ------- + ∇ × f : numpy.ndarray + Curl of f. + """ - n = _np.prod(field.shape[3:]) - k_s = _ks(size,field.shape[:3],True) + n = _np.prod(f.shape[3:]) + k_s = _ks(size,f.shape[:3],True) e = _np.zeros((3, 3, 3)) e[0, 1, 2] = e[1, 2, 0] = e[2, 0, 1] = +1.0 # Levi-Civita symbol e[0, 2, 1] = e[2, 1, 0] = e[1, 0, 2] = -1.0 - field_fourier = _np.fft.rfftn(field,axes=(0,1,2)) - curl_ = (_np.einsum('slm,ijkl,ijkm ->ijks', e,k_s,field_fourier)*2.0j*_np.pi if n == 3 else # vector, 3 -> 3 - _np.einsum('slm,ijkl,ijknm->ijksn',e,k_s,field_fourier)*2.0j*_np.pi) # tensor, 3x3 -> 3x3 + f_fourier = _np.fft.rfftn(f,axes=(0,1,2)) + curl_ = (_np.einsum('slm,ijkl,ijkm ->ijks', e,k_s,f_fourier)*2.0j*_np.pi if n == 3 else # vector, 3 -> 3 + _np.einsum('slm,ijkl,ijknm->ijksn',e,k_s,f_fourier)*2.0j*_np.pi) # tensor, 3x3 -> 3x3 - return _np.fft.irfftn(curl_,axes=(0,1,2),s=field.shape[:3]) + return _np.fft.irfftn(curl_,axes=(0,1,2),s=f.shape[:3]) -def divergence(size,field): - """ +def divergence(size,f): + u""" Calculate divergence of a vector or tensor field in Fourier space. Parameters ---------- size : numpy.ndarray of shape (3) Physical size of the periodic field. - field : numpy.ndarray of shape (:,:,:,3) or (:,:,:,3,3) + f : numpy.ndarray of shape (:,:,:,3) or (:,:,:,3,3) Periodic field of which the divergence is calculated. + Returns + ------- + ∇ · f : numpy.ndarray + Divergence of f. + """ - n = _np.prod(field.shape[3:]) - k_s = _ks(size,field.shape[:3],True) + n = _np.prod(f.shape[3:]) + k_s = _ks(size,f.shape[:3],True) - field_fourier = _np.fft.rfftn(field,axes=(0,1,2)) - div_ = (_np.einsum('ijkl,ijkl ->ijk', k_s,field_fourier)*2.0j*_np.pi if n == 3 else # vector, 3 -> 1 - _np.einsum('ijkm,ijklm->ijkl',k_s,field_fourier)*2.0j*_np.pi) # tensor, 3x3 -> 3 + f_fourier = _np.fft.rfftn(f,axes=(0,1,2)) + div_ = (_np.einsum('ijkl,ijkl ->ijk', k_s,f_fourier)*2.0j*_np.pi if n == 3 else # vector, 3 -> 1 + _np.einsum('ijkm,ijklm->ijkl',k_s,f_fourier)*2.0j*_np.pi) # tensor, 3x3 -> 3 - return _np.fft.irfftn(div_,axes=(0,1,2),s=field.shape[:3]) + return _np.fft.irfftn(div_,axes=(0,1,2),s=f.shape[:3]) -def gradient(size,field): - """ - Calculate gradient of a scalar or vector field in Fourier space. +def gradient(size,f): + u""" + Calculate gradient of a scalar or vector fieldin Fourier space. Parameters ---------- size : numpy.ndarray of shape (3) Physical size of the periodic field. - field : numpy.ndarray of shape (:,:,:,1) or (:,:,:,3) + f : numpy.ndarray of shape (:,:,:,1) or (:,:,:,3) Periodic field of which the gradient is calculated. + Returns + ------- + ∇ f : numpy.ndarray + Divergence of f. + """ - n = _np.prod(field.shape[3:]) - k_s = _ks(size,field.shape[:3],True) + n = _np.prod(f.shape[3:]) + k_s = _ks(size,f.shape[:3],True) - field_fourier = _np.fft.rfftn(field,axes=(0,1,2)) - grad_ = (_np.einsum('ijkl,ijkm->ijkm', field_fourier,k_s)*2.0j*_np.pi if n == 1 else # scalar, 1 -> 3 - _np.einsum('ijkl,ijkm->ijklm',field_fourier,k_s)*2.0j*_np.pi) # vector, 3 -> 3x3 + f_fourier = _np.fft.rfftn(f,axes=(0,1,2)) + grad_ = (_np.einsum('ijkl,ijkm->ijkm', f_fourier,k_s)*2.0j*_np.pi if n == 1 else # scalar, 1 -> 3 + _np.einsum('ijkl,ijkm->ijklm',f_fourier,k_s)*2.0j*_np.pi) # vector, 3 -> 3x3 - return _np.fft.irfftn(grad_,axes=(0,1,2),s=field.shape[:3]) + return _np.fft.irfftn(grad_,axes=(0,1,2),s=f.shape[:3]) def coordinates0_point(cells,size,origin=_np.zeros(3)): @@ -123,6 +139,11 @@ def coordinates0_point(cells,size,origin=_np.zeros(3)): origin : numpy.ndarray, optional Physical origin of the periodic field. Defaults to [0.0,0.0,0.0]. + Returns + ------- + x_p_0 : numpy.ndarray + Undeformed cell center coordinates. + """ start = origin + size/cells*.5 end = origin + size - size/cells*.5 @@ -144,6 +165,11 @@ def displacement_fluct_point(size,F): F : numpy.ndarray Deformation gradient field. + Returns + ------- + u_p_fluct : numpy.ndarray + Fluctuating part of the cell center displacements. + """ integrator = 0.5j*size/_np.pi @@ -171,6 +197,11 @@ def displacement_avg_point(size,F): F : numpy.ndarray Deformation gradient field. + Returns + ------- + u_p_avg : numpy.ndarray + Average part of the cell center displacements. + """ F_avg = _np.average(F,axis=(0,1,2)) return _np.einsum('ml,ijkl->ijkm',F_avg - _np.eye(3),coordinates0_point(F.shape[:3],size)) @@ -187,6 +218,11 @@ def displacement_point(size,F): F : numpy.ndarray Deformation gradient field. + Returns + ------- + u_p : numpy.ndarray + Cell center displacements. + """ return displacement_avg_point(size,F) + displacement_fluct_point(size,F) @@ -204,6 +240,11 @@ def coordinates_point(size,F,origin=_np.zeros(3)): origin : numpy.ndarray of shape (3), optional Physical origin of the periodic field. Defaults to [0.0,0.0,0.0]. + Returns + ------- + x_p : numpy.ndarray + Cell center coordinates. + """ return coordinates0_point(F.shape[:3],size,origin) + displacement_point(size,F) @@ -220,6 +261,11 @@ def cellsSizeOrigin_coordinates0_point(coordinates0,ordered=True): Expect coordinates0 data to be ordered (x fast, z slow). Defaults to True. + Returns + ------- + cells, size, origin : Three numpy.ndarray, each of shape (3) + Information to reconstruct grid. + """ coords = [_np.unique(coordinates0[:,i]) for i in range(3)] mincorner = _np.array(list(map(min,coords))) @@ -252,19 +298,6 @@ def cellsSizeOrigin_coordinates0_point(coordinates0,ordered=True): return (cells,size,origin) -def coordinates0_check(coordinates0): - """ - Check whether coordinates lie on a regular grid. - - Parameters - ---------- - coordinates0 : numpy.ndarray - Array of undeformed cell coordinates. - - """ - cellsSizeOrigin_coordinates0_point(coordinates0,ordered=True) - - def coordinates0_node(cells,size,origin=_np.zeros(3)): """ Nodal positions (undeformed). @@ -278,6 +311,11 @@ def coordinates0_node(cells,size,origin=_np.zeros(3)): origin : numpy.ndarray of shape (3), optional Physical origin of the periodic field. Defaults to [0.0,0.0,0.0]. + Returns + ------- + x_n_0 : numpy.ndarray + Undeformed nodal coordinates. + """ return _np.stack(_np.meshgrid(_np.linspace(origin[0],size[0]+origin[0],cells[0]+1), _np.linspace(origin[1],size[1]+origin[1],cells[1]+1), @@ -296,6 +334,11 @@ def displacement_fluct_node(size,F): F : numpy.ndarray Deformation gradient field. + Returns + ------- + u_n_fluct : numpy.ndarray + Fluctuating part of the nodal displacements. + """ return point_to_node(displacement_fluct_point(size,F)) @@ -311,6 +354,11 @@ def displacement_avg_node(size,F): F : numpy.ndarray Deformation gradient field. + Returns + ------- + u_n_avg : numpy.ndarray + Average part of the nodal displacements. + """ F_avg = _np.average(F,axis=(0,1,2)) return _np.einsum('ml,ijkl->ijkm',F_avg - _np.eye(3),coordinates0_node(F.shape[:3],size)) @@ -327,6 +375,11 @@ def displacement_node(size,F): F : numpy.ndarray Deformation gradient field. + Returns + ------- + u_p : numpy.ndarray + Nodal displacements. + """ return displacement_avg_node(size,F) + displacement_fluct_node(size,F) @@ -344,28 +397,15 @@ def coordinates_node(size,F,origin=_np.zeros(3)): origin : numpy.ndarray of shape (3), optional Physical origin of the periodic field. Defaults to [0.0,0.0,0.0]. + Returns + ------- + x_n : numpy.ndarray + Nodal coordinates. + """ return coordinates0_node(F.shape[:3],size,origin) + displacement_node(size,F) -def point_to_node(cell_data): - """Interpolate periodic point data to nodal data.""" - n = ( cell_data + _np.roll(cell_data,1,(0,1,2)) - + _np.roll(cell_data,1,(0,)) + _np.roll(cell_data,1,(1,)) + _np.roll(cell_data,1,(2,)) - + _np.roll(cell_data,1,(0,1)) + _np.roll(cell_data,1,(1,2)) + _np.roll(cell_data,1,(2,0)))*0.125 - - return _np.pad(n,((0,1),(0,1),(0,1))+((0,0),)*len(cell_data.shape[3:]),mode='wrap') - - -def node_2_point(node_data): - """Interpolate periodic nodal data to point data.""" - c = ( node_data + _np.roll(node_data,1,(0,1,2)) - + _np.roll(node_data,1,(0,)) + _np.roll(node_data,1,(1,)) + _np.roll(node_data,1,(2,)) - + _np.roll(node_data,1,(0,1)) + _np.roll(node_data,1,(1,2)) + _np.roll(node_data,1,(2,0)))*0.125 - - return c[1:,1:,1:] - - def cellsSizeOrigin_coordinates0_node(coordinates0,ordered=True): """ Return grid 'DNA', i.e. cells, size, and origin from 1D array of nodal positions. @@ -378,6 +418,11 @@ def cellsSizeOrigin_coordinates0_node(coordinates0,ordered=True): Expect coordinates0 data to be ordered (x fast, z slow). Defaults to True. + Returns + ------- + cells, size, origin : Three numpy.ndarray, each of shape (3) + Information to reconstruct grid. + """ coords = [_np.unique(coordinates0[:,i]) for i in range(3)] mincorner = _np.array(list(map(min,coords))) @@ -402,6 +447,72 @@ def cellsSizeOrigin_coordinates0_node(coordinates0,ordered=True): return (cells,size,origin) +def point_to_node(cell_data): + """ + Interpolate periodic point data to nodal data. + + Parameters + ---------- + cell_data : numpy.ndarray of shape (:,:,:,...) + Data defined on the cell centers of a periodic grid. + + Returns + ------- + node_data : numpy.ndarray of shape (:,:,:,...) + Data defined on the nodes of a periodic grid. + + """ + n = ( cell_data + _np.roll(cell_data,1,(0,1,2)) + + _np.roll(cell_data,1,(0,)) + _np.roll(cell_data,1,(1,)) + _np.roll(cell_data,1,(2,)) + + _np.roll(cell_data,1,(0,1)) + _np.roll(cell_data,1,(1,2)) + _np.roll(cell_data,1,(2,0)))*0.125 + + return _np.pad(n,((0,1),(0,1),(0,1))+((0,0),)*len(cell_data.shape[3:]),mode='wrap') + + +def node_to_point(node_data): + """ + Interpolate periodic nodal data to point data. + + Parameters + ---------- + node_data : numpy.ndarray of shape (:,:,:,...) + Data defined on the nodes of a periodic grid. + + Returns + ------- + cell_data : numpy.ndarray of shape (:,:,:,...) + Data defined on the cell centers of a periodic grid. + + """ + c = ( node_data + _np.roll(node_data,1,(0,1,2)) + + _np.roll(node_data,1,(0,)) + _np.roll(node_data,1,(1,)) + _np.roll(node_data,1,(2,)) + + _np.roll(node_data,1,(0,1)) + _np.roll(node_data,1,(1,2)) + _np.roll(node_data,1,(2,0)))*0.125 + + return c[1:,1:,1:] + + +def coordinates0_valid(coordinates0): + """ + Check whether coordinates lie on a regular grid. + + Parameters + ---------- + coordinates0 : numpy.ndarray + Array of undeformed cell coordinates. + + Returns + ------- + valid : bool + Wheter the coordinates lie on a regular grid. + + """ + try: + cellsSizeOrigin_coordinates0_point(coordinates0,ordered=True) + return True + except ValueError: + return False + + def regrid(size,F,cells): """ Return mapping from coordinates in deformed configuration to a regular grid. diff --git a/python/damask/mechanics.py b/python/damask/mechanics.py index 66f5c9916..6ba785df1 100644 --- a/python/damask/mechanics.py +++ b/python/damask/mechanics.py @@ -1,4 +1,9 @@ -"""Finite-strain continuum mechanics.""" +""" +Finite-strain continuum mechanics. + +All routines operate on numpy.ndarrays of shape (...,3,3). + +""" from . import tensor as _tensor from . import _rotation @@ -154,7 +159,6 @@ def strain(F,t,m): return eps - def stress_Cauchy(P,F): """ Calculate the Cauchy stress (true stress). diff --git a/python/damask/seeds.py b/python/damask/seeds.py index 9ab148a82..aeab2e64a 100644 --- a/python/damask/seeds.py +++ b/python/damask/seeds.py @@ -24,6 +24,11 @@ def from_random(size,N_seeds,cells=None,rng_seed=None): A seed to initialize the BitGenerator. Defaults to None. If None, then fresh, unpredictable entropy will be pulled from the OS. + Returns + ------- + coords : numpy.ndarray of shape (N_seeds,3) + Seed coordinates in 3D space. + """ rng = _np.random.default_rng(rng_seed) if cells is None: @@ -56,6 +61,11 @@ def from_Poisson_disc(size,N_seeds,N_candidates,distance,periodic=True,rng_seed= A seed to initialize the BitGenerator. Defaults to None. If None, then fresh, unpredictable entropy will be pulled from the OS. + Returns + ------- + coords : numpy.ndarray of shape (N_seeds,3) + Seed coordinates in 3D space. + """ rng = _np.random.default_rng(rng_seed) coords = _np.empty((N_seeds,3)) @@ -94,6 +104,11 @@ def from_grid(grid,selection=None,invert=False,average=False,periodic=True): periodic : boolean, optional Center of gravity with periodic boundaries. + Returns + ------- + coords, materials : numpy.ndarray of shape (:,3), numpy.ndarray of shape (:) + Seed coordinates in 3D space, material IDs. + """ material = grid.material.reshape((-1,1),order='F') mask = _np.full(grid.cells.prod(),True,dtype=bool) if selection is None else \ diff --git a/python/damask/solver/__init__.py b/python/damask/solver/__init__.py index 298adb1f7..01cfeed9b 100644 --- a/python/damask/solver/__init__.py +++ b/python/damask/solver/__init__.py @@ -1,3 +1,3 @@ -"""Tools to control the various solvers.""" +"""Run simulations directly from python.""" from ._marc import Marc # noqa diff --git a/python/damask/tensor.py b/python/damask/tensor.py index 00fc7e72a..cf5d94020 100644 --- a/python/damask/tensor.py +++ b/python/damask/tensor.py @@ -1,10 +1,7 @@ """ -Tensor operations. +Tensor mathematics. -Notes ------ -This is not a tensor class, but a collection of routines -to operate on numpy.ndarrays of shape (...,3,3). +All routines operate on numpy.ndarrays of shape (...,3,3). """ diff --git a/python/damask/util.py b/python/damask/util.py index b519fb15e..fcf41c251 100644 --- a/python/damask/util.py +++ b/python/damask/util.py @@ -49,7 +49,7 @@ _colors = { #################################################################################################### def srepr(arg,glue = '\n'): r""" - Join arguments with glue string. + Join items with glue string. Parameters ---------- @@ -58,6 +58,11 @@ def srepr(arg,glue = '\n'): glue : str, optional Glue used for joining operation. Defaults to \n. + Returns + ------- + joined : str + String representation of the joined items. + """ if (not hasattr(arg, 'strip') and (hasattr(arg, '__getitem__') or @@ -68,19 +73,71 @@ def srepr(arg,glue = '\n'): def emph(what): - """Formats string with emphasis.""" + """ + Format with emphasis. + + Parameters + ---------- + what : object with __repr__ or iterable of objects with __repr__. + Message to format. + + Returns + ------- + formatted : str + Formatted string representation of the joined items. + + """ return _colors['bold']+srepr(what)+_colors['end_color'] def deemph(what): - """Formats string with deemphasis.""" + """ + Format with deemphasis. + + Parameters + ---------- + what : object with __repr__ or iterable of objects with __repr__. + Message to format. + + Returns + ------- + formatted : str + Formatted string representation of the joined items. + + """ return _colors['dim']+srepr(what)+_colors['end_color'] def warn(what): - """Formats string for warning.""" + """ + Format for warning. + + Parameters + ---------- + what : object with __repr__ or iterable of objects with __repr__. + Message to format. + + Returns + ------- + formatted : str + Formatted string representation of the joined items. + + """ return _colors['warning']+emph(what)+_colors['end_color'] def strikeout(what): - """Formats string as strikeout.""" + """ + Format as strikeout. + + Parameters + ---------- + what : object with __repr__ or iterable of objects with __repr__. + Message to format. + + Returns + ------- + formatted : str + Formatted string representation of the joined items. + + """ return _colors['crossout']+srepr(what)+_colors['end_color'] @@ -97,6 +154,11 @@ def execute(cmd,wd='./',env=None): env : dict, optional Environment for execution. + Returns + ------- + stdout, stderr : str + Output of the executed command. + """ print(f"executing '{cmd}' in '{wd}'") process = subprocess.run(shlex.split(cmd), @@ -157,6 +219,11 @@ def scale_to_coprime(v): v : numpy.ndarray of shape (:) Vector to scale. + Returns + ------- + m : numpy.ndarray of shape (:) + Vector scaled to co-prime numbers. + """ MAX_DENOMINATOR = 1000000 @@ -371,9 +438,14 @@ def DREAM3D_base_group(fname): Parameters ---------- - fname : str + fname : str or pathlib.Path Filename of the DREAM.3D (HDF5) file. + Returns + ------- + path : str + Path to the base group. + """ with h5py.File(fname,'r') as f: base_group = f.visit(lambda path: path.rsplit('/',2)[0] if '_SIMPL_GEOMETRY/SPACING' in path else None) @@ -393,9 +465,14 @@ def DREAM3D_cell_data_group(fname): Parameters ---------- - fname : str + fname : str or pathlib.Path Filename of the DREAM.3D (HDF5) file. + Returns + ------- + path : str + Path to the cell data group. + """ base_group = DREAM3D_base_group(fname) with h5py.File(fname,'r') as f: diff --git a/python/tests/test_Result.py b/python/tests/test_Result.py index 0c34e7989..e79e25784 100644 --- a/python/tests/test_Result.py +++ b/python/tests/test_Result.py @@ -271,7 +271,7 @@ class TestResult: @pytest.mark.parametrize('overwrite',['off','on']) def test_add_overwrite(self,default,overwrite): - last = default.view('times',default.times_in_range(0,np.inf)[-1]) + last = default.view('increments',-1) last.add_stress_Cauchy() @@ -279,9 +279,9 @@ class TestResult: created_first = datetime.strptime(created_first,'%Y-%m-%d %H:%M:%S%z') if overwrite == 'on': - last = last.allow_modification() + last = last.modification_enable() else: - last = last.disallow_modification() + last = last.modification_disable() time.sleep(2.) try: @@ -301,14 +301,24 @@ class TestResult: def test_rename(self,default,allowed): if allowed == 'on': F = default.place('F') - default = default.allow_modification() + default = default.modification_enable() default.rename('F','new_name') assert np.all(F == default.place('new_name')) - default = default.disallow_modification() + default = default.modification_disable() with pytest.raises(PermissionError): default.rename('P','another_new_name') + @pytest.mark.parametrize('allowed',['off','on']) + def test_remove(self,default,allowed): + if allowed == 'on': + unsafe = default.modification_enable() + unsafe.remove('F') + assert unsafe.get('F') is None + else: + with pytest.raises(PermissionError): + default.remove('F') + @pytest.mark.parametrize('mode',['cell','node']) def test_coordinates(self,default,mode): if mode == 'cell': diff --git a/python/tests/test_grid_filters.py b/python/tests/test_grid_filters.py index d43e94c3c..d9c074c11 100644 --- a/python/tests/test_grid_filters.py +++ b/python/tests/test_grid_filters.py @@ -60,7 +60,7 @@ class TestGridFilters: cell_field_x = np.interp(coordinates0_point_x,coordinates_node_x,node_field_x,period=np.pi*2.) cell_field = np.broadcast_to(cell_field_x.reshape(-1,1,1),cells) - assert np.allclose(cell_field,grid_filters.node_2_point(node_field)) + assert np.allclose(cell_field,grid_filters.node_to_point(node_field)) @pytest.mark.parametrize('mode',['point','node']) def test_coordinates0_origin(self,mode): @@ -93,14 +93,24 @@ class TestGridFilters: F = np.broadcast_to(np.random.random((3,3)), tuple(cells)+(3,3)) assert np.allclose(function(size,F),0.0) - @pytest.mark.parametrize('function',[grid_filters.coordinates0_check, - grid_filters.cellsSizeOrigin_coordinates0_node, - grid_filters.cellsSizeOrigin_coordinates0_point]) + @pytest.mark.parametrize('function',[grid_filters.cellsSizeOrigin_coordinates0_point, + grid_filters.cellsSizeOrigin_coordinates0_node]) def test_invalid_coordinates(self,function): invalid_coordinates = np.random.random((np.random.randint(12,52),3)) with pytest.raises(ValueError): function(invalid_coordinates) + @pytest.mark.parametrize('function',[grid_filters.coordinates0_point, + grid_filters.coordinates0_node]) + def test_valid_coordinates_check(self,function): + valid_coordinates = function(np.random.randint(4,10,(3)),np.random.rand(3)) + assert grid_filters.coordinates0_valid(valid_coordinates.reshape(-1,3,order='F')) + + def test_invalid_coordinates_check(self): + invalid_coordinates = np.random.random((np.random.randint(12,52),3)) + assert not grid_filters.coordinates0_valid(invalid_coordinates) + + @pytest.mark.parametrize('function',[grid_filters.cellsSizeOrigin_coordinates0_node, grid_filters.cellsSizeOrigin_coordinates0_point]) def test_uneven_spaced_coordinates(self,function): diff --git a/src/commercialFEM_fileList.f90 b/src/commercialFEM_fileList.f90 index 1eac4f8fc..5ddc7219e 100644 --- a/src/commercialFEM_fileList.f90 +++ b/src/commercialFEM_fileList.f90 @@ -40,7 +40,6 @@ #include "phase_damage_isobrittle.f90" #include "phase_damage_isoductile.f90" #include "phase_damage_anisobrittle.f90" -#include "phase_damage_anisoductile.f90" #include "homogenization.f90" #include "homogenization_mechanical.f90" #include "homogenization_mechanical_pass.f90" diff --git a/src/material.f90 b/src/material.f90 index 5802ab3c8..add9fd1c3 100644 --- a/src/material.f90 +++ b/src/material.f90 @@ -35,7 +35,7 @@ module material integer, dimension(:,:), allocatable, public, protected :: & ! (constituent,elem) material_phaseAt, & !< phase ID of each element material_phaseID, & !< per (constituent,cell) - material_phaseEntry !< per (constituent,cell + material_phaseEntry !< per (constituent,cell) integer, dimension(:,:,:), allocatable, public, protected :: & ! (constituent,IP,elem) material_phaseMemberAt !< position of the element within its phase instance diff --git a/src/phase_damage.f90 b/src/phase_damage.f90 index 6d2243933..5a1d5970b 100644 --- a/src/phase_damage.f90 +++ b/src/phase_damage.f90 @@ -12,8 +12,7 @@ submodule(phase) damage DAMAGE_UNDEFINED_ID, & DAMAGE_ISOBRITTLE_ID, & DAMAGE_ISODUCTILE_ID, & - DAMAGE_ANISOBRITTLE_ID, & - DAMAGE_ANISODUCTILE_ID + DAMAGE_ANISOBRITTLE_ID end enum @@ -34,10 +33,6 @@ submodule(phase) damage logical, dimension(:), allocatable :: mySources end function anisobrittle_init - module function anisoductile_init() result(mySources) - logical, dimension(:), allocatable :: mySources - end function anisoductile_init - module function isobrittle_init() result(mySources) logical, dimension(:), allocatable :: mySources end function isobrittle_init @@ -62,10 +57,6 @@ submodule(phase) damage S end subroutine anisobrittle_dotState - module subroutine anisoductile_dotState(ph,me) - integer, intent(in) :: ph,me - end subroutine anisoductile_dotState - module subroutine isoductile_dotState(ph,me) integer, intent(in) :: ph,me end subroutine isoductile_dotState @@ -75,11 +66,6 @@ submodule(phase) damage character(len=*), intent(in) :: group end subroutine anisobrittle_results - module subroutine anisoductile_results(phase,group) - integer, intent(in) :: phase - character(len=*), intent(in) :: group - end subroutine anisoductile_results - module subroutine isobrittle_results(phase,group) integer, intent(in) :: phase character(len=*), intent(in) :: group @@ -146,7 +132,6 @@ module subroutine damage_init where(isobrittle_init() ) phase_source = DAMAGE_ISOBRITTLE_ID where(isoductile_init() ) phase_source = DAMAGE_ISODUCTILE_ID where(anisobrittle_init()) phase_source = DAMAGE_ANISOBRITTLE_ID - where(anisoductile_init()) phase_source = DAMAGE_ANISODUCTILE_ID endif end subroutine damage_init @@ -171,7 +156,7 @@ module function phase_f_phi(phi,co,ce) result(f) en = material_phaseEntry(co,ce) select case(phase_source(ph)) - case(DAMAGE_ISOBRITTLE_ID,DAMAGE_ISODUCTILE_ID,DAMAGE_ANISOBRITTLE_ID,DAMAGE_ANISODUCTILE_ID) + case(DAMAGE_ISOBRITTLE_ID,DAMAGE_ISODUCTILE_ID,DAMAGE_ANISOBRITTLE_ID) f = 1.0_pReal & - phi*damageState(ph)%state(1,en) case default @@ -304,9 +289,6 @@ module subroutine damage_results(group,ph) case (DAMAGE_ANISOBRITTLE_ID) sourceType call anisobrittle_results(ph,group//'damage/') - case (DAMAGE_ANISODUCTILE_ID) sourceType - call anisoductile_results(ph,group//'damage/') - end select sourceType end subroutine damage_results @@ -332,9 +314,6 @@ function phase_damage_collectDotState(ph,me) result(broken) case (DAMAGE_ISODUCTILE_ID) sourceType call isoductile_dotState(ph,me) - case (DAMAGE_ANISODUCTILE_ID) sourceType - call anisoductile_dotState(ph,me) - case (DAMAGE_ANISOBRITTLE_ID) sourceType call anisobrittle_dotState(mechanical_S(ph,me), ph,me) ! correct stress? diff --git a/src/phase_damage_anisoductile.f90 b/src/phase_damage_anisoductile.f90 deleted file mode 100644 index 049f148bc..000000000 --- a/src/phase_damage_anisoductile.f90 +++ /dev/null @@ -1,138 +0,0 @@ -!-------------------------------------------------------------------------------------------------- -!> @author Luv Sharma, Max-Planck-Institut für Eisenforschung GmbH -!> @author Pratheek Shanthraj, Max-Planck-Institut für Eisenforschung GmbH -!> @brief material subroutine incorporating anisotropic ductile damage source mechanism -!> @details to be done -!-------------------------------------------------------------------------------------------------- -submodule(phase:damage) anisoductile - - type :: tParameters !< container type for internal constitutive parameters - real(pReal) :: & - q !< damage rate sensitivity - real(pReal), dimension(:), allocatable :: & - gamma_crit !< critical plastic strain per slip system - character(len=pStringLen), allocatable, dimension(:) :: & - output - end type tParameters - - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters - -contains - - -!-------------------------------------------------------------------------------------------------- -!> @brief module initialization -!> @details reads in material parameters, allocates arrays, and does sanity checks -!-------------------------------------------------------------------------------------------------- -module function anisoductile_init() result(mySources) - - logical, dimension(:), allocatable :: mySources - - class(tNode), pointer :: & - phases, & - phase, & - mech, & - pl, & - sources, & - src - integer :: Ninstances,Nmembers,ph - integer, dimension(:), allocatable :: N_sl - character(len=pStringLen) :: extmsg = '' - - - mySources = source_active('anisoductile') - if(count(mySources) == 0) return - - print'(/,a)', ' <<<+- phase:damage:anisoductile init -+>>>' - print'(a,i0)', ' # phases: ',count(mySources); flush(IO_STDOUT) - - - phases => config_material%get('phase') - allocate(param(phases%length)) - - do ph = 1, phases%length - if(mySources(ph)) then - phase => phases%get(ph) - mech => phase%get('mechanical') - pl => mech%get('plastic') - sources => phase%get('damage') - - - associate(prm => param(ph)) - src => sources%get(1) - - N_sl = pl%get_as1dInt('N_sl',defaultVal=emptyIntArray) - prm%q = src%get_asFloat('q') - prm%gamma_crit = src%get_as1dFloat('gamma_crit',requiredSize=size(N_sl)) - - ! expand: family => system - prm%gamma_crit = math_expand(prm%gamma_crit,N_sl) - -#if defined (__GFORTRAN__) - prm%output = output_as1dString(src) -#else - prm%output = src%get_as1dString('output',defaultVal=emptyStringArray) -#endif - - ! sanity checks - if (prm%q <= 0.0_pReal) extmsg = trim(extmsg)//' q' - if (any(prm%gamma_crit < 0.0_pReal)) extmsg = trim(extmsg)//' gamma_crit' - - Nmembers=count(material_phaseID==ph) - call phase_allocateState(damageState(ph),Nmembers,1,1,0) - damageState(ph)%atol = src%get_asFloat('anisoDuctile_atol',defaultVal=1.0e-3_pReal) - if(any(damageState(ph)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' anisoductile_atol' - - end associate - -!-------------------------------------------------------------------------------------------------- -! exit if any parameter is out of range - if (extmsg /= '') call IO_error(211,ext_msg=trim(extmsg)//'(damage_anisoDuctile)') - endif - - enddo - - -end function anisoductile_init - - -!-------------------------------------------------------------------------------------------------- -!> @brief calculates derived quantities from state -!-------------------------------------------------------------------------------------------------- -module subroutine anisoductile_dotState(ph,me) - - integer, intent(in) :: & - ph, & - me - - - associate(prm => param(ph)) - damageState(ph)%dotState(1,me) = sum(plasticState(ph)%slipRate(:,me)/(damage_phi(ph,me)**prm%q)/prm%gamma_crit) - end associate - -end subroutine anisoductile_dotState - - -!-------------------------------------------------------------------------------------------------- -!> @brief writes results to HDF5 output file -!-------------------------------------------------------------------------------------------------- -module subroutine anisoductile_results(phase,group) - - integer, intent(in) :: phase - character(len=*), intent(in) :: group - - integer :: o - - - associate(prm => param(phase), stt => damageState(phase)%state) - outputsLoop: do o = 1,size(prm%output) - select case(trim(prm%output(o))) - case ('f_phi') - call results_writeDataset(group,stt,trim(prm%output(o)),'driving force','J/m³') - end select - enddo outputsLoop - end associate - -end subroutine anisoductile_results - -end submodule anisoductile diff --git a/src/phase_thermal.f90 b/src/phase_thermal.f90 index f6d046266..b51ac74f5 100644 --- a/src/phase_thermal.f90 +++ b/src/phase_thermal.f90 @@ -190,7 +190,7 @@ end function phase_thermal_collectDotState !-------------------------------------------------------------------------------------------------- -!> @brief Damage viscosity. +!> @brief Thermal viscosity. !-------------------------------------------------------------------------------------------------- module function phase_mu_T(co,ce) result(mu) @@ -205,7 +205,7 @@ end function phase_mu_T !-------------------------------------------------------------------------------------------------- -!> @brief Damage conductivity/diffusivity in reference configuration. +!> @brief Thermal conductivity/diffusivity in reference configuration. !-------------------------------------------------------------------------------------------------- module function phase_K_T(co,ce) result(K)