Merge branch 'development' of magit1.mpie.de:damask/DAMASK into development

DAMASK_env.sh

@@ -2,7 +2,7 @@
 # usage:  source DAMASK_env.sh
 
 if [ "$OSTYPE" == "linux-gnu" ] || [ "$OSTYPE" == 'linux' ]; then
-  DAMASK_ROOT=$(readlink -f "`dirname $BASH_SOURCE`")
+  DAMASK_ROOT=$(python -c "import os,sys; print(os.path.realpath(os.path.expanduser(sys.argv[1])))" "`dirname $BASH_SOURCE`")
 else
   [[ "${BASH_SOURCE::1}" == "/" ]] && BASE="" || BASE="`pwd`/"
   STAT=$(stat "`dirname $BASE$BASH_SOURCE`")
@@ -18,11 +18,11 @@ fi
 
 SOLVER=`which DAMASK_spectral 2>/dev/null`
 if [ "x$SOLVER" == "x" ]; then
-  export SOLVER='Not found!'
+  SOLVER='Not found!'
 fi
 PROCESSING=`which postResults 2>/dev/null`
 if [ "x$PROCESSING" == "x" ]; then
-  export PROCESSING='Not found!'
+  PROCESSING='Not found!'
 fi
 
 # according to http://software.intel.com/en-us/forums/topic/501500
@@ -55,7 +55,8 @@ if [ ! -z "$PS1" ]; then
   echo "Multithreading     DAMASK_NUM_THREADS=$DAMASK_NUM_THREADS"
   if [ "x$PETSC_DIR"   != "x" ]; then
     echo "PETSc location     $PETSC_DIR"
-    [[ `readlink -f $PETSC_DIR` == $PETSC_DIR ]] || echo "               ~~> "`readlink -f $PETSC_DIR`
+    [[ `python -c "import os,sys; print(os.path.realpath(os.path.expanduser(sys.argv[1])))" "$PETSC_DIR"` == $PETSC_DIR ]] \
+    || echo "               ~~> "`python -c "import os,sys; print(os.path.realpath(os.path.expanduser(sys.argv[1])))" "$PETSC_DIR"`
   fi
   [[ "x$PETSC_ARCH"  != "x" ]] && echo "PETSc architecture $PETSC_ARCH"
   echo "MSC.Marc/Mentat    $MSC_ROOT"

code/DAMASK_spectral.f90

@@ -13,7 +13,8 @@ program DAMASK_spectral
    pInt, &
    pLongInt, &
    pReal, &
-   tol_math_check
+   tol_math_check, &
+   dNeq
  use DAMASK_interface, only: &
    DAMASK_interface_init, &
    loadCaseFile, &
@@ -147,7 +148,9 @@ program DAMASK_spectral
    MPI_file_seek, &
    MPI_file_get_position, &
    MPI_file_write, &
-   MPI_allreduce
+   MPI_abort, &
+   MPI_allreduce, &
+   PETScFinalize
 
 !--------------------------------------------------------------------------------------------------
 ! init DAMASK (all modules)
@@ -339,7 +342,7 @@ program DAMASK_spectral
                 reshape(spread(tol_math_check,1,9),[ 3,3]))&
                 .or. abs(math_det33(loadCases(currentLoadCase)%rotation)) > &
                 1.0_pReal + tol_math_check) errorID = 846_pInt                                      ! given rotation matrix contains strain
-     if (any(loadCases(currentLoadCase)%rotation /= math_I3)) &
+     if (any(dNeq(loadCases(currentLoadCase)%rotation, math_I3))) &
        write(6,'(2x,a,/,3(3(3x,f12.7,1x)/))',advance='no') 'rotation of loadframe:',&
                 math_transpose33(loadCases(currentLoadCase)%rotation)
      if (loadCases(currentLoadCase)%time < 0.0_pReal)          errorID = 834_pInt                   ! negative time increment
@@ -423,17 +426,21 @@ program DAMASK_spectral
 !--------------------------------------------------------------------------------------------------
 ! prepare MPI parallel out (including opening of file)
  allocate(outputSize(worldsize), source = 0_MPI_OFFSET_KIND)
- outputSize(worldrank+1) = int(size(materialpoint_results)*pReal,MPI_OFFSET_KIND)
+ outputSize(worldrank+1) = size(materialpoint_results,kind=MPI_OFFSET_KIND)*int(pReal,MPI_OFFSET_KIND)
  call MPI_allreduce(MPI_IN_PLACE,outputSize,worldsize,MPI_INT,MPI_SUM,PETSC_COMM_WORLD,ierr)        ! get total output size over each process
+ if(ierr /=0_pInt) call IO_error(894_pInt, ext_msg='MPI_allreduce')
  call MPI_file_open(PETSC_COMM_WORLD, &
                     trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//'.spectralOut', &
                     MPI_MODE_WRONLY + MPI_MODE_APPEND, &
                     MPI_INFO_NULL, &
                     resUnit, &
                     ierr)
+ if(ierr /=0_pInt) call IO_error(894_pInt, ext_msg='MPI_file_open')
  call MPI_file_get_position(resUnit,fileOffset,ierr)                                                ! get offset from header
+ if(ierr /=0_pInt) call IO_error(894_pInt, ext_msg='MPI_file_get_position')
  fileOffset = fileOffset + sum(outputSize(1:worldrank))                                             ! offset of my process in file (header + processes before me)
  call MPI_file_seek (resUnit,fileOffset,MPI_SEEK_SET,ierr)
+ if(ierr /=0_pInt) call IO_error(894_pInt, ext_msg='MPI_file_seek')
 
  if (.not. appendToOutFile) then                                                                    ! if not restarting, write 0th increment
    do i=1, size(materialpoint_results,3)/(maxByteOut/(materialpoint_sizeResults*pReal))+1           ! slice the output of my process in chunks not exceeding the limit for one output
@@ -443,6 +450,7 @@ program DAMASK_spectral
                                    [(outputIndex(2)-outputIndex(1)+1)*materialpoint_sizeResults]), &
                          (outputIndex(2)-outputIndex(1)+1)*materialpoint_sizeResults,&
                          MPI_DOUBLE, MPI_STATUS_IGNORE, ierr)
+     if(ierr /=0_pInt) call IO_error(894_pInt, ext_msg='MPI_file_write')
    enddo
    fileOffset = fileOffset + sum(outputSize)                                                        ! forward to current file position
    if (worldrank == 0) &
@@ -643,6 +651,7 @@ program DAMASK_spectral
            write(6,'(1/,a)') ' ... writing results to file ......................................'
          call materialpoint_postResults()
          call MPI_file_seek (resUnit,fileOffset,MPI_SEEK_SET,ierr)
+         if(ierr /=0_pInt) call IO_error(894_pInt, ext_msg='MPI_file_seek')
          do i=1, size(materialpoint_results,3)/(maxByteOut/(materialpoint_sizeResults*pReal))+1     ! slice the output of my process in chunks not exceeding the limit for one output
            outputIndex=int([(i-1_pInt)*((maxByteOut/pReal)/materialpoint_sizeResults)+1_pInt, &
                       min(i*((maxByteOut/pReal)/materialpoint_sizeResults),size(materialpoint_results,3))],pLongInt)
@@ -650,6 +659,7 @@ program DAMASK_spectral
                                          [(outputIndex(2)-outputIndex(1)+1)*materialpoint_sizeResults]), &
                                (outputIndex(2)-outputIndex(1)+1)*materialpoint_sizeResults,&
                                MPI_DOUBLE, MPI_STATUS_IGNORE, ierr)
+           if(ierr /=0_pInt) call IO_error(894_pInt, ext_msg='MPI_file_write')
          enddo
          fileOffset = fileOffset + sum(outputSize)                                                  ! forward to current file position
        endif
@@ -698,7 +708,7 @@ program DAMASK_spectral
  enddo
  call utilities_destroy()
 
- call PetscFinalize(ierr); CHKERRQ(ierr)
+ call PETScFinalize(ierr); CHKERRQ(ierr)
 
  if (notConvergedCounter > 0_pInt) call quit(3_pInt)                                                ! error if some are not converged
  call quit(0_pInt)                                                                                  ! no complains ;)

code/IO.f90

@@ -1669,6 +1669,8 @@ subroutine IO_error(error_ID,el,ip,g,ext_msg)
    msg = 'unknown filter type selected'
  case (893_pInt)
    msg = 'PETSc: SNES_DIVERGED_FNORM_NAN'
+ case (894_pInt)
+   msg = 'MPI error'
  
 !-------------------------------------------------------------------------------------------------
 ! error messages related to parsing of Abaqus input file

code/material.f90

@@ -1280,7 +1280,7 @@ subroutine material_populateGrains
  integer(pInt) :: t,e,i,g,j,m,c,r,homog,micro,sgn,hme, myDebug, &
                   phaseID,textureID,dGrains,myNgrains,myNorientations,myNconstituents, &
                   grain,constituentGrain,ipGrain,symExtension, ip
- real(pReal) :: extreme,rnd
+ real(pReal) :: deviation,extreme,rnd
  integer(pInt),  dimension (:,:),   allocatable :: Nelems                                           ! counts number of elements in homog, micro array
  type(p_intvec), dimension (:,:), allocatable :: elemsOfHomogMicro                                  ! lists element number in homog, micro array
 
@@ -1407,8 +1407,11 @@ subroutine material_populateGrains
          extreme = 0.0_pReal
          t = 0_pInt
          do i = 1_pInt,myNconstituents                                                              ! find largest deviator
-           if (real(sgn,pReal)*log(NgrainsOfConstituent(i)/myNgrains/microstructure_fraction(i,micro)) > extreme) then
+           deviation = real(sgn,pReal)*log( microstructure_fraction(i,micro) / &
-             extreme = real(sgn,pReal)*log(NgrainsOfConstituent(i)/myNgrains/microstructure_fraction(i,micro))
+                                           !-------------------------------- &
+                                           (real(NgrainsOfConstituent(i),pReal)/real(myNgrains,pReal) ) )
+           if (deviation > extreme) then
+             extreme = deviation
              t = i
            endif
          enddo

code/prec.f90

@@ -113,7 +113,9 @@ module prec
 
  public :: &
    prec_init, &
-   prec_isNaN
+   prec_isNaN, &
+   dEq, &
+   dNeq
  
 contains
 

config/Phase_DisloKMC_Tungsten.config

@@ -1,37 +0,0 @@
-### $Id$ ###
-
-[Tungsten]
-elasticity          hooke
-plasticity          dislokmc
-
-### Material parameters ###
-lattice_structure     bcc
-C11                   523.0e9       # From Marinica et al. Journal of Physics: Condensed Matter(2013)
-C12                   202.0e9
-C44                   161.0e9
-
-grainsize             2.0e-5        # Average grain size [m] 2.0e-5
-SolidSolutionStrength 0.0           # Strength due to elements in solid solution
-
-### Dislocation glide parameters ###
-#per family
-Nslip               12 0 
-slipburgers         2.72e-10        # Burgers vector of slip system [m]
-rhoedge0            1.0e12          # Initial edge dislocation density [m/m**3]
-rhoedgedip0         1.0             # Initial edged dipole dislocation density [m/m**3]
-Qedge               2.725e-19       # Activation energy for dislocation glide [J]
-v0                  3560.3          # Initial glide velocity [m/s](kmC)
-p_slip              0.16            # p-exponent in glide velocity
-q_slip              1.00            # q-exponent in glide velocity                                                                  
-u_slip              2.47            # u-exponent of stress pre-factor (kmC)                                                                     
-s_slip              0.97            # self hardening in glide velocity (kmC)
-tau_peierls         2.03e9          # peierls stress [Pa]
-
-#hardening
-dipoleformationfactor 0             # to have hardening due to dipole formation off
-CLambdaSlip         10.0            # Adj. parameter controlling dislocation mean free path
-D0                  4.0e-5          # Vacancy diffusion prefactor [m**2/s]
-Qsd                 4.5e-19         # Activation energy for climb [J]
-Catomicvolume       1.0             # Adj. parameter controlling the atomic volume [in b]
-Cedgedipmindistance 1.0             # Adj. parameter controlling the minimum dipole distance [in b]
-interaction_slipslip	 0.2 0.11 0.19 0.15 0.11 0.17

lib/damask/asciitable.py

@@ -40,12 +40,12 @@ class ASCIItable():
       self.__IO__['in'] = name
 
     try:
-      self.__IO__['out'] = (open(outname,'w') if (not os.path.isfile(outname) \
+      self.__IO__['out'] = (open(outname,'w') if (not os.path.isfile(outname) or
-                                                   or os.access(     outname, os.W_OK) \
+                                                      os.access(     outname, os.W_OK)
-                                                 ) \
+                                                 ) and
-                                              and (not self.__IO__['inPlace'] \
+                                                 (not self.__IO__['inPlace'] or
-                                                   or not os.path.isfile(name) \
+                                                  not os.path.isfile(name)   or
-                                                   or     os.access(     name, os.W_OK) \
+                                                      os.access(     name, os.W_OK)
                                                  ) else None) if outname else sys.stdout
     except TypeError:
       self.__IO__['out'] = outname

processing/post/averageTable.py

@@ -37,9 +37,13 @@ if options.label is None:
 if filenames == []: filenames = [None]
 
 for name in filenames:
-  try:
+  damask.util.croak(name)
-    table = damask.ASCIItable(name = name,
+
-                              outname = options.label+'_averaged_'+name if name else name,
+  try:    table = damask.ASCIItable(name = name,
+                                    outname = os.path.join(
+                                              os.path.split(name)[0],
+                                              options.label+'_averaged_'+os.path.split(name)[1]
+                                              ) if name else name,
                                     buffered = False)
   except: continue
   damask.util.report(scriptName,name)

processing/post/histogram.py

@@ -34,7 +34,11 @@ parser.add_option('-N',
                   dest = 'N',
                   type = 'int', metavar = 'int',
                   help = 'number of bins')
-parser.add_option('-l', '--logarithmic',
+parser.add_option('--density',
+                  dest = 'density',
+                  action = 'store_true',
+                  help = 'report probability density')
+parser.add_option('--logarithmic',
                   dest = 'log',
                   action = 'store_true',
                   help = 'logarithmically spaced bins')
@@ -42,6 +46,7 @@ parser.set_defaults(data    = None,
                     weights = None,
                     range   = None,
                     N       = None,
+                    density = False,
                     log     = False,
                    )
 
@@ -110,6 +115,7 @@ for name in filenames:
                                            (forward(rangeMax)-forward(rangeMin))/options.N),
                              range = (rangeMin,rangeMax),
                              weights = None if options.weights is None else table.data[:,1],
+                             density = options.density,
                             )
   bincenter = reverse(forward(rangeMin) + (0.5+np.arange(options.N)) *
                      (forward(rangeMax)-forward(rangeMin))/options.N)                             # determine center of bins

processing/pre/geom_fromVoronoiTessellation.py

@@ -4,7 +4,7 @@
 import os,sys,math
 import numpy as np
 import multiprocessing
-from optparse import OptionParser
+from optparse import OptionParser,OptionGroup
 from scipy    import spatial
 import damask
 
@@ -109,71 +109,83 @@ Generate geometry description and material configuration by standard Voronoi tes
 
 """, version = scriptID)
 
-parser.add_option('-g', '--grid',
+
-                  dest = 'grid',
+group = OptionGroup(parser, "Tessellation","")
-                  type = 'int', nargs = 3, metavar = ' '.join(['int']*3),
+
-                  help = 'a,b,c grid of hexahedral box [auto]')
+group.add_option('-l', '--laguerre',
-parser.add_option('-s', '--size',
-                  dest = 'size',
-                  type = 'float', nargs = 3, metavar=' '.join(['float']*3),
-                  help = 'x,y,z size of hexahedral box [auto]')
-parser.add_option('-o', '--origin',
-                  dest = 'origin',
-                  type = 'float', nargs = 3, metavar=' '.join(['float']*3),
-                  help = 'offset from old to new origin of grid')
-parser.add_option('-p', '--position',
-                  dest = 'position',
-                  type = 'string', metavar = 'string',
-                  help = 'column label for seed positions [%default]')
-parser.add_option('-w', '--weight',
-                  dest = 'weight',
-                  type = 'string', metavar = 'string',
-                  help = 'column label for seed weights [%default]')
-parser.add_option('-m', '--microstructure',
-                  dest = 'microstructure',
-                  type = 'string', metavar = 'string',
-                  help = 'column label for seed microstructures [%default]')
-parser.add_option('-e', '--eulers',
-                  dest = 'eulers',
-                  type = 'string', metavar = 'string',
-                  help = 'column label for seed Euler angles [%default]')
-parser.add_option('--axes',
-                  dest = 'axes',
-                  type = 'string', nargs = 3, metavar = ' '.join(['string']*3),
-                  help = 'orientation coordinate frame in terms of position coordinate frame')
-parser.add_option('--homogenization',
-                  dest = 'homogenization',
-                  type = 'int', metavar = 'int',
-                  help = 'homogenization index to be used [%default]')
-parser.add_option('--crystallite',
-                  dest = 'crystallite',
-                  type = 'int', metavar = 'int',
-                  help = 'crystallite index to be used [%default]')
-parser.add_option('--phase',
-                  dest = 'phase',
-                  type = 'int', metavar = 'int',
-                  help = 'phase index to be used [%default]')
-parser.add_option('-r', '--rnd',
-                  dest = 'randomSeed',
-                  type = 'int', metavar='int',
-                  help = 'seed of random number generator for second phase distribution [%default]')
-parser.add_option('--secondphase',
-                  dest = 'secondphase',
-                  type = 'float', metavar= 'float',
-                  help = 'volume fraction of randomly distribute second phase [%default]')
-parser.add_option('-l', '--laguerre',
                   dest = 'laguerre',
                   action = 'store_true',
                   help = 'use Laguerre (weighted Voronoi) tessellation')
-parser.add_option('--cpus',
+group.add_option('--cpus',
                   dest = 'cpus',
                   type = 'int', metavar = 'int',
                   help = 'number of parallel processes to use for Laguerre tessellation [%default]')
-parser.add_option('--nonperiodic',
+group.add_option('--nonperiodic',
                   dest = 'nonperiodic',
                   action = 'store_true',
                   help = 'use nonperiodic tessellation')
 
+parser.add_option_group(group)
+
+group = OptionGroup(parser, "Geometry","")
+
+group.add_option('-g', '--grid',
+                  dest = 'grid',
+                  type = 'int', nargs = 3, metavar = ' '.join(['int']*3),
+                  help = 'a,b,c grid of hexahedral box [auto]')
+group.add_option('-s', '--size',
+                  dest = 'size',
+                  type = 'float', nargs = 3, metavar=' '.join(['float']*3),
+                  help = 'x,y,z size of hexahedral box [auto]')
+group.add_option('-o', '--origin',
+                  dest = 'origin',
+                  type = 'float', nargs = 3, metavar=' '.join(['float']*3),
+                  help = 'origin of grid')
+
+parser.add_option_group(group)
+
+group = OptionGroup(parser, "Seeds","")
+
+group.add_option('-p', '--position',
+                  dest = 'position',
+                  type = 'string', metavar = 'string',
+                  help = 'column label for seed positions [%default]')
+group.add_option('-w', '--weight',
+                  dest = 'weight',
+                  type = 'string', metavar = 'string',
+                  help = 'column label for seed weights [%default]')
+group.add_option('-m', '--microstructure',
+                  dest = 'microstructure',
+                  type = 'string', metavar = 'string',
+                  help = 'column label for seed microstructures [%default]')
+group.add_option('-e', '--eulers',
+                  dest = 'eulers',
+                  type = 'string', metavar = 'string',
+                  help = 'column label for seed Euler angles [%default]')
+group.add_option('--axes',
+                  dest = 'axes',
+                  type = 'string', nargs = 3, metavar = ' '.join(['string']*3),
+                  help = 'orientation coordinate frame in terms of position coordinate frame')
+
+parser.add_option_group(group)
+
+group = OptionGroup(parser, "Configuration","")
+
+group.add_option('--homogenization',
+                  dest = 'homogenization',
+                  type = 'int', metavar = 'int',
+                  help = 'homogenization index to be used [%default]')
+group.add_option('--crystallite',
+                  dest = 'crystallite',
+                  type = 'int', metavar = 'int',
+                  help = 'crystallite index to be used [%default]')
+group.add_option('--phase',
+                  dest = 'phase',
+                  type = 'int', metavar = 'int',
+                  help = 'phase index to be used [%default]')
+
+parser.add_option_group(group)
+
 parser.set_defaults(position       = 'pos',
                     weight         = 'weight',
                     microstructure = 'microstructure',
@@ -181,24 +193,18 @@ parser.set_defaults(position       = 'pos',
                     homogenization = 1,
                     crystallite    = 1,
                     phase          = 1,
-                    secondphase    = 0.0,
                     cpus           = 2,
                     laguerre       = False,
                     nonperiodic    = False,
-                    randomSeed     = None,
                   )
 (options,filenames) = parser.parse_args()
 
-if options.secondphase > 1.0 or options.secondphase < 0.0:
- parser.error('volume fraction of second phase ({}) out of bounds.'.format(options.secondphase))
-
 # --- loop over input files -------------------------------------------------------------------------
 
 if filenames == []: filenames = [None]
 
 for name in filenames:
-  try:
+  try:    table = damask.ASCIItable(name = name,
-    table = damask.ASCIItable(name = name,
                                     outname = os.path.splitext(name)[-2]+'.geom' if name else name,
                                     buffered = False)
   except: continue
@@ -294,20 +300,11 @@ for name in filenames:
   config_header = []
   formatwidth = 1+int(math.log10(NgrainIDs))
 
-  phase = options.phase * np.ones(NgrainIDs,'i')
-  if int(options.secondphase*info['microstructures']) > 0:
-    phase[0:int(options.secondphase*info['microstructures'])] += 1                                # alter fraction 'options.secondphase' of used grainIDs
-    randomSeed = options.randomSeed if options.randomSeed \
-                 else int(os.urandom(4).encode('hex'), 16)                                        # random seed for second phase
-    np.random.seed(randomSeed)
-    np.random.shuffle(phase)
-    config_header += ['# random seed (phase shuffling): {}'.format(randomSeed)]
-
   config_header += ['<microstructure>']
   for i,ID in enumerate(grainIDs):
     config_header += ['[Grain%s]'%(str(ID).zfill(formatwidth)),
                       'crystallite %i'%options.crystallite,
-                      '(constituent)\tphase %i\ttexture %s\tfraction 1.0'%(phase[i],str(ID).rjust(formatwidth)),
+                      '(constituent)\tphase %i\ttexture %s\tfraction 1.0'%(options.phase,str(ID).rjust(formatwidth)),
                      ]
   if hasEulers:
     config_header += ['<texture>']

processing/pre/seeds_fromRandom.py

@@ -48,8 +48,11 @@ parser.add_option('-m', '--microstructure',
 parser.add_option('-r', '--rnd',
                   dest = 'randomSeed', type = 'int', metavar = 'int',
                   help = 'seed of random number generator [%default]')
+parser.add_option('--format',
+                  dest = 'format', type = 'string', metavar = 'string',
+                  help = 'number format of output [auto]')
 
-group = OptionGroup(parser, "Laguerre Tessellation Options",
+group = OptionGroup(parser, "Laguerre Tessellation",
                    "Parameters determining shape of weight distribution of seed points"
                    )
 group.add_option('-w', '--weights',
@@ -70,8 +73,8 @@ group.add_option('--sigma',
                  help='standard deviation of normal distribution for weights [%default]')
 parser.add_option_group(group)
 
-group = OptionGroup(parser, "Selective Seeding Options",
+group = OptionGroup(parser, "Selective Seeding",
-                    "More uniform distribution of seed points using Mitchell\'s Best Candidate Algorithm"
+                    "More uniform distribution of seed points using Mitchell's Best Candidate Algorithm"
                    )
 group.add_option('-s','--selective',
                   action = 'store_true',
@@ -103,6 +106,7 @@ parser.set_defaults(randomSeed = None,
                     force = False,
                     distance = 0.2,
                     numCandidates = 10,
+                    format = None,
                    )
 
 (options,filenames) = parser.parse_args()
@@ -215,7 +219,7 @@ for name in filenames:
 # --- write seeds information ------------------------------------------------------------
 
   table.data = seeds
-  table.data_writeArray()
+  table.data_writeArray(fmt = options.format)
     
 # --- output finalization --------------------------------------------------------------------------