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

This commit is contained in:
Philip Eisenlohr 2016-09-27 14:24:24 -04:00
commit 03b023a0d1
13 changed files with 177 additions and 236 deletions

View File

@ -1 +1 @@
v2.0.1-148-gfb23c49
v2.0.1-150-g5345b42

View File

@ -9,6 +9,14 @@
/* http://stackoverflow.com/questions/30279228/is-there-an-alternative-to-getcwd-in-fortran-2003-2008 */
int isdirectory_c(const char *dir){
struct stat statbuf;
if(stat(dir, &statbuf) != 0)
return 0;
return S_ISDIR(statbuf.st_mode);
}
void getcurrentworkdir_c(char cwd[], int *stat ){
char cwd_tmp[1024];
if(getcwd(cwd_tmp, sizeof(cwd_tmp)) == cwd_tmp){
@ -20,9 +28,14 @@ void getcurrentworkdir_c(char cwd[], int *stat ){
}
}
int isdirectory_c(const char *dir){
struct stat statbuf;
if(stat(dir, &statbuf) != 0)
return 0;
return S_ISDIR(statbuf.st_mode);
void gethostname_c(char hostname[], int *stat ){
char hostname_tmp[1024];
if(gethostname(hostname_tmp, sizeof(hostname_tmp)) == 0){
strcpy(hostname,hostname_tmp);
*stat = 0;
}
else{
*stat = 1;
}
}

View File

@ -1568,8 +1568,6 @@ subroutine IO_error(error_ID,el,ip,g,ext_msg)
msg = 'math_check: R*v == q*v failed'
case (410_pInt)
msg = 'eigenvalues computation error'
case (450_pInt)
msg = 'unknown symmetry type specified'
!-------------------------------------------------------------------------------------------------
! homogenization errors

View File

@ -1864,37 +1864,29 @@ end function math_sampleGaussVar
!--------------------------------------------------------------------------------------------------
!> @brief symmetrically equivalent Euler angles for given sample symmetry 1:triclinic, 2:monoclinic, 4:orthotropic
!> @brief symmetrically equivalent Euler angles for given sample symmetry
!> @detail 1 (equivalent to != 2 and !=4):triclinic, 2:monoclinic, 4:orthotropic
!--------------------------------------------------------------------------------------------------
pure function math_symmetricEulers(sym,Euler)
implicit none
integer(pInt), intent(in) :: sym
integer(pInt), intent(in) :: sym !< symmetry Class
real(pReal), dimension(3), intent(in) :: Euler
real(pReal), dimension(3,3) :: math_symmetricEulers
integer(pInt) :: i,j
math_symmetricEulers(1,1) = PI+Euler(1)
math_symmetricEulers(2,1) = Euler(2)
math_symmetricEulers(3,1) = Euler(3)
math_symmetricEulers = transpose(reshape([PI+Euler(1), PI-Euler(1), 2.0_pReal*PI-Euler(1), &
Euler(2), PI-Euler(2), PI -Euler(2), &
Euler(3), PI+Euler(3), PI +Euler(3)],[3,3])) ! transpose is needed to have symbolic notation instead of column-major
math_symmetricEulers(1,2) = PI-Euler(1)
math_symmetricEulers(2,2) = PI-Euler(2)
math_symmetricEulers(3,2) = PI+Euler(3)
math_symmetricEulers(1,3) = 2.0_pReal*PI-Euler(1)
math_symmetricEulers(2,3) = PI-Euler(2)
math_symmetricEulers(3,3) = PI+Euler(3)
forall (i=1_pInt:3_pInt,j=1_pInt:3_pInt) math_symmetricEulers(j,i) = modulo(math_symmetricEulers(j,i),2.0_pReal*pi)
math_symmetricEulers = modulo(math_symmetricEulers,2.0_pReal*pi)
select case (sym)
case (4_pInt) ! all done
case (4_pInt) ! orthotropic: all done
case (2_pInt) ! return only first
case (2_pInt) ! monoclinic: return only first
math_symmetricEulers(1:3,2:3) = 0.0_pReal
case default ! return blank
case default ! triclinic: return blank
math_symmetricEulers = 0.0_pReal
end select

View File

@ -20,6 +20,7 @@ module DAMASK_interface
geometryFile = '', & !< parameter given for geometry file
loadCaseFile = '' !< parameter given for load case file
character(len=1024), private :: workingDirectory !< accessed by getSolverWorkingDirectoryName for compatibility reasons
character, private,parameter :: pathSep = '/'
public :: &
getSolverWorkingDirectoryName, &
@ -31,7 +32,6 @@ module DAMASK_interface
getLoadCaseFile, &
rectifyPath, &
makeRelativePath, &
getPathSep, &
IIO_stringValue, &
IIO_intValue, &
IIO_lc, &
@ -44,6 +44,8 @@ contains
!--------------------------------------------------------------------------------------------------
subroutine DAMASK_interface_init()
use, intrinsic :: iso_fortran_env ! to get compiler_version and compiler_options (at least for gfortran 4.6 at the moment)
use system_routines, only: &
getHostName
implicit none
character(len=1024) :: &
@ -64,6 +66,7 @@ subroutine DAMASK_interface_init()
integer, dimension(8) :: &
dateAndTime ! type default integer
PetscErrorCode :: ierr
logical :: error
external :: &
quit,&
MPI_Comm_rank,&
@ -116,7 +119,6 @@ subroutine DAMASK_interface_init()
tag = IIO_lc(IIO_stringValue(commandLine,chunkPos,i)) ! extract key
select case(tag)
case ('-h','--help')
mainProcess2: if (worldrank == 0) then
write(6,'(a)') ' #######################################################################'
write(6,'(a)') ' DAMASK_spectral:'
write(6,'(a)') ' The spectral method boundary value problem solver for'
@ -163,7 +165,6 @@ subroutine DAMASK_interface_init()
write(6,'(/,a)')' --help'
write(6,'(a,/)')' Prints this message and exits'
call quit(0_pInt) ! normal Termination
endif mainProcess2
case ('-l', '--load', '--loadcase')
loadcaseArg = IIO_stringValue(commandLine,chunkPos,i+1_pInt)
case ('-g', '--geom', '--geometry')
@ -185,12 +186,11 @@ subroutine DAMASK_interface_init()
geometryFile = getGeometryFile(geometryArg)
loadCaseFile = getLoadCaseFile(loadCaseArg)
call get_environment_variable('HOSTNAME',hostName)
call get_environment_variable('USER',userName)
mainProcess3: if (worldrank == 0) then
error = getHostName(hostName)
write(6,'(a,a)') ' Host name: ', trim(hostName)
write(6,'(a,a)') ' User name: ', trim(userName)
write(6,'(a,a)') ' Path separator: ', getPathSep()
write(6,'(a,a)') ' Path separator: ', pathSep
write(6,'(a,a)') ' Command line call: ', trim(commandLine)
if (len(trim(workingDirArg))>0) &
write(6,'(a,a)') ' Working dir argument: ', trim(workingDirArg)
@ -203,7 +203,6 @@ subroutine DAMASK_interface_init()
if (SpectralRestartInc > 1_pInt) &
write(6,'(a,i6.6)') ' Restart at increment: ', spectralRestartInc
write(6,'(a,l1,/)') ' Append to result file: ', appendToOutFile
endif mainProcess3
end subroutine DAMASK_interface_init
@ -222,11 +221,9 @@ character(len=1024) function storeWorkingDirectory(workingDirectoryArg,geometryA
character(len=*), intent(in) :: workingDirectoryArg !< working directory argument
character(len=*), intent(in) :: geometryArg !< geometry argument
character(len=1024) :: cwd
character :: pathSep
logical :: error
external :: quit
pathSep = getPathSep()
wdGiven: if (len(workingDirectoryArg)>0) then
absolutePath: if (workingDirectoryArg(1:1) == pathSep) then
storeWorkingDirectory = workingDirectoryArg
@ -262,6 +259,7 @@ end function storeWorkingDirectory
character(len=1024) function getSolverWorkingDirectoryName()
implicit none
getSolverWorkingDirectoryName = workingDirectory
end function getSolverWorkingDirectoryName
@ -274,10 +272,8 @@ character(len=1024) function getSolverJobName()
implicit none
integer :: posExt,posSep
character :: pathSep
character(len=1024) :: tempString
pathSep = getPathSep()
tempString = geometryFile
posExt = scan(tempString,'.',back=.true.)
@ -308,11 +304,9 @@ character(len=1024) function getGeometryFile(geometryParameter)
cwd
integer :: posExt, posSep
logical :: error
character :: pathSep
external :: quit
getGeometryFile = geometryParameter
pathSep = getPathSep()
posExt = scan(getGeometryFile,'.',back=.true.)
posSep = scan(getGeometryFile,pathSep,back=.true.)
@ -344,11 +338,9 @@ character(len=1024) function getLoadCaseFile(loadCaseParameter)
cwd
integer :: posExt, posSep
logical :: error
character :: pathSep
external :: quit
getLoadCaseFile = loadcaseParameter
pathSep = getPathSep()
posExt = scan(getLoadCaseFile,'.',back=.true.)
posSep = scan(getLoadCaseFile,pathSep,back=.true.)
@ -374,11 +366,8 @@ function rectifyPath(path)
implicit none
character(len=*) :: path
character(len=len_trim(path)) :: rectifyPath
character :: pathSep
integer :: i,j,k,l ! no pInt
pathSep = getPathSep()
!--------------------------------------------------------------------------------------------------
! remove /./ from path
l = len_trim(path)
@ -415,10 +404,8 @@ character(len=1024) function makeRelativePath(a,b)
implicit none
character (len=*) :: a,b
character :: pathSep
integer :: i,posLastCommonSlash,remainingSlashes !no pInt
pathSep = getPathSep()
posLastCommonSlash = 0
remainingSlashes = 0
@ -434,35 +421,6 @@ character(len=1024) function makeRelativePath(a,b)
end function makeRelativePath
!--------------------------------------------------------------------------------------------------
!> @brief counting / and \ in $PATH System variable the character occuring more often is assumed
! to be the path separator
!--------------------------------------------------------------------------------------------------
character function getPathSep()
implicit none
character(len=2048) :: &
path
integer(pInt) :: &
backslash = 0_pInt, &
slash = 0_pInt
integer :: i
call get_environment_variable('PATH',path)
do i=1, len(trim(path))
if (path(i:i)=='/') slash = slash + 1_pInt
if (path(i:i)=='\') backslash = backslash + 1_pInt
enddo
if (backslash>slash) then
getPathSep = '\'
else
getPathSep = '/'
endif
end function getPathSep
!--------------------------------------------------------------------------------------------------
!> @brief taken from IO, check IO_stringValue for documentation
!--------------------------------------------------------------------------------------------------

View File

@ -9,7 +9,8 @@ module system_routines
public :: &
isDirectory, &
getCWD
getCWD, &
getHostName
interface
@ -29,6 +30,14 @@ interface
integer(C_INT),intent(out) :: stat
end subroutine getCurrentWorkDir_C
subroutine getHostName_C(str, stat) bind(C)
use, intrinsic :: ISO_C_Binding, only: &
C_INT, &
C_CHAR
character(kind=C_CHAR), dimension(1024), intent(out) :: str ! C string is an array
integer(C_INT),intent(out) :: stat
end subroutine getHostName_C
end interface
@ -85,5 +94,34 @@ logical function getCWD(str)
end function getCWD
!--------------------------------------------------------------------------------------------------
!> @brief gets the current host name
!--------------------------------------------------------------------------------------------------
logical function getHostName(str)
use, intrinsic :: ISO_C_Binding, only: &
C_INT, &
C_CHAR, &
C_NULL_CHAR
implicit none
character(len=*), intent(out) :: str
character(kind=C_CHAR), dimension(1024) :: strFixedLength ! C string is an array
integer(C_INT) :: stat
integer :: i
str = repeat('',len(str))
call getHostName_C(strFixedLength,stat)
do i=1,1024 ! copy array components until Null string is found
if (strFixedLength(i) /= C_NULL_CHAR) then
str(i:i)=strFixedLength(i)
else
exit
endif
enddo
getHostName=merge(.True.,.False.,stat /= 0_C_INT)
end function getHostName
end module system_routines

View File

@ -8,32 +8,20 @@ plasticity nonlocal
(output) rho_edge
(output) rho_screw
(output) rho_sgl
(output) rho_sgl_edge
(output) rho_sgl_edge_pos
(output) rho_sgl_edge_neg
(output) rho_sgl_screw
(output) rho_sgl_screw_pos
(output) rho_sgl_screw_neg
(output) rho_sgl_mobile
(output) rho_sgl_edge_mobile
(output) rho_sgl_edge_pos_mobile
(output) rho_sgl_edge_neg_mobile
(output) rho_sgl_screw_mobile
(output) rho_sgl_screw_pos_mobile
(output) rho_sgl_screw_neg_mobile
(output) rho_sgl_immobile
(output) rho_sgl_edge_immobile
(output) rho_sgl_edge_pos_immobile
(output) rho_sgl_edge_neg_immobile
(output) rho_sgl_screw_immobile
(output) rho_sgl_screw_pos_immobile
(output) rho_sgl_screw_neg_immobile
(output) rho_dip
(output) rho_dip_edge
(output) rho_dip_screw
(output) excess_rho
(output) excess_rho_edge
(output) excess_rho_screw
(output) rho_forest
(output) delta
(output) delta_sgl
@ -47,10 +35,8 @@ plasticity nonlocal
(output) rho_dot_sgl
(output) rho_dot_sgl_mobile
(output) rho_dot_dip
(output) rho_dot_gen
(output) rho_dot_gen_edge
(output) rho_dot_gen_screw
(output) rho_dot_sgl2dip
(output) rho_dot_sgl2dip_edge
(output) rho_dot_sgl2dip_screw
(output) rho_dot_ann_ath
@ -66,24 +52,6 @@ plasticity nonlocal
(output) velocity_edge_neg
(output) velocity_screw_pos
(output) velocity_screw_neg
(output) slipdirection.x
(output) slipdirection.y
(output) slipdirection.z
(output) slipnormal.x
(output) slipnormal.y
(output) slipnormal.z
(output) fluxDensity_edge_pos.x
(output) fluxDensity_edge_pos.y
(output) fluxDensity_edge_pos.z
(output) fluxDensity_edge_neg.x
(output) fluxDensity_edge_neg.y
(output) fluxDensity_edge_neg.z
(output) fluxDensity_screw_pos.x
(output) fluxDensity_screw_pos.y
(output) fluxDensity_screw_pos.z
(output) fluxDensity_screw_neg.x
(output) fluxDensity_screw_neg.y
(output) fluxDensity_screw_neg.z
(output) maximumDipoleHeight_edge
(output) maximumDipoleHeight_screw
(output) accumulated_shear

20
examples/ConfigFiles/Phase_Nonlocal_Nickel.config Executable file → Normal file
View File

@ -14,8 +14,6 @@ plasticity nonlocal
(output) rho_dip_edge
(output) rho_dip_screw
(output) rho_forest
(output) excess_rho_edge
(output) excess_rho_screw
(output) accumulatedshear
(output) shearrate
(output) resolvedstress
@ -30,24 +28,6 @@ plasticity nonlocal
(output) rho_dot_edgejogs
(output) rho_dot_flux_edge
(output) rho_dot_flux_screw
(output) slipdirection.x
(output) slipdirection.y
(output) slipdirection.z
(output) slipnormal.x
(output) slipnormal.y
(output) slipnormal.z
(output) fluxdensity_edge_pos.x
(output) fluxdensity_edge_pos.y
(output) fluxdensity_edge_pos.z
(output) fluxdensity_edge_neg.x
(output) fluxdensity_edge_neg.y
(output) fluxdensity_edge_neg.z
(output) fluxdensity_screw_pos.x
(output) fluxdensity_screw_pos.y
(output) fluxdensity_screw_pos.z
(output) fluxdensity_screw_neg.x
(output) fluxdensity_screw_neg.y
(output) fluxdensity_screw_neg.z
lattice_structure fcc
Nslip 12 # number of slip systems per family

View File

@ -7,11 +7,6 @@ plasticity phenopowerlaw
(output) resolvedstress_slip
(output) accumulated_shear_slip
(output) totalshear
(output) resistance_twin
(output) shearrate_twin
(output) resolvedstress_twin
(output) accumulated_shear_twin
(output) totalvolfrac_twin
lattice_structure fcc
Nslip 12 # per family
@ -26,19 +21,6 @@ n_slip 20
tau0_slip 31e6 # per family
tausat_slip 63e6 # per family
a_slip 2.25
gdot0_twin 0.001
n_twin 20
tau0_twin 31e6 # per family
s_pr 0 # push-up factor for slip saturation due to twinning
twin_b 0
twin_c 0
twin_d 0
twin_e 0
h0_slipslip 75e6
h0_twinslip 0
h0_twintwin 0
interaction_slipslip 1 1 1.4 1.4 1.4 1.4
interaction_sliptwin 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
interaction_twinslip 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
interaction_twintwin 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
atol_resistance 1

View File

@ -23,6 +23,8 @@ class Test():
'keep': False,
'accept': False,
'updateRequest': False,
'show': False,
'select': None,
}
for arg in defaults.keys():
setattr(self,arg,kwargs.get(arg) if kwargs.get(arg) else defaults[arg])
@ -58,10 +60,18 @@ class Test():
action = "store_true",
dest = "accept",
help = "calculate results but always consider test as successfull")
self.parser.add_option("-l", "--list",
action = "store_true",
dest = "show",
help = "show all test variants and do no calculation")
self.parser.add_option("-s", "--select",
dest = "select",
help = "run test of given name only")
self.parser.set_defaults(keep = self.keep,
accept = self.accept,
update = self.updateRequest,
show = self.show,
select = self.select,
)
@ -73,6 +83,11 @@ class Test():
self.prepareAll()
for variant,name in enumerate(self.variants):
if self.options.show:
logging.critical('{}: {}'.format(variant,name))
elif self.options.select is not None and name != self.options.select:
pass
else:
try:
if not self.options.keep:
self.prepare(variant)

View File

@ -11,7 +11,7 @@ do
< $geom \
| \
vtk_addRectilinearGridData \
--scalar microstructure \
--data microstructure \
--inplace \
--vtk ${geom%.*}.vtk
rm ${geom%.*}.vtk

View File

@ -270,7 +270,6 @@ def rcbParser(content,M,size,tolerance,idcolumn,segmentcolumn):
myWalk = point['segments'].pop()
grainLegs = [myWalk]
myEnd = segments[myWalk][1 if segments[myWalk][0] == myStart else 0]
while (myEnd != pointId):
myV = [points[myEnd]['coords'][0]-points[myStart]['coords'][0],
points[myEnd]['coords'][1]-points[myStart]['coords'][1]]
@ -281,13 +280,13 @@ def rcbParser(content,M,size,tolerance,idcolumn,segmentcolumn):
if peek == myWalk:
continue # do not go back same path
peekEnd = segments[peek][1 if segments[peek][0] == myEnd else 0]
peekV = [points[myEnd]['coords'][0]-points[peekEnd]['coords'][0],
points[myEnd]['coords'][1]-points[peekEnd]['coords'][1]]
peekV = [points[peekEnd]['coords'][0]-points[myEnd]['coords'][0],
points[peekEnd]['coords'][1]-points[myEnd]['coords'][1]]
peekLen = math.sqrt(peekV[0]**2+peekV[1]**2)
if peekLen == 0.0: damask.util.croak('peeklen is zero: peek point {}'.format(peek))
crossproduct = (myV[0]*peekV[1] - myV[1]*peekV[0])/myLen/peekLen
dotproduct = (myV[0]*peekV[0] + myV[1]*peekV[1])/myLen/peekLen
innerAngle = crossproduct*(dotproduct+1.0)
innerAngle = math.copysign(1.0,crossproduct)*(dotproduct-1.0)
if innerAngle >= best['product']: # takes sharpest left turn
best['product'] = innerAngle
best['peek'] = peek
@ -301,7 +300,7 @@ def rcbParser(content,M,size,tolerance,idcolumn,segmentcolumn):
if myWalk in points[myStart]['segments']:
points[myStart]['segments'].remove(myWalk)
else:
damask.utilcroak('{} not in segments of point {}'.format(myWalk,myStart))
damask.util.croak('{} not in segments of point {}'.format(myWalk,myStart))
grainDraw.append(points[myStart]['coords'])
grainLegs.append(myWalk)
@ -311,7 +310,6 @@ def rcbParser(content,M,size,tolerance,idcolumn,segmentcolumn):
else:
grains['box'] = grainLegs
# build overall data structure
rcData = {'dimension':[dX,dY],
@ -772,8 +770,8 @@ def fftbuild(rcData,height,xframe,yframe,resolution,extrusion):
'dimension':(xsize,ysize,zsize)}
frameindex=len(rcData['grain'])+1 # calculate frame index as largest grain index plus one
dx = xsize/(xres+1) # calculate step sizes
dy = ysize/(yres+1)
dx = xsize/(xres) # calculate step sizes
dy = ysize/(yres)
grainpoints = []
for segments in rcData['grain']: # get segments of each grain
@ -788,11 +786,10 @@ def fftbuild(rcData,height,xframe,yframe,resolution,extrusion):
grainpoints.append([]) # start out blank for current grain
for p in sorted(points, key=points.get): # loop thru set of sorted points
grainpoints[-1].append([rcData['point'][p][0],rcData['point'][p][1]]) # append x,y of point
bestGuess = 0 # assume grain 0 as best guess
for i in range(int(xres*yres)): # walk through all points in xy plane
xtest = -xframe+((i%xres)+0.5)*dx # calculate coordinates
ytest = -yframe+(int(i/xres)+0.5)*dy
ytest = -yframe+((i//xres)+0.5)*dy
if(xtest < 0 or xtest > maxX): # check wether part of frame
if( ytest < 0 or ytest > maxY): # part of edges
fftdata['fftpoints'].append(frameindex+2) # append frameindex to result array

View File

@ -5,7 +5,7 @@ do
vtk_pointcloud $seeds
vtk_addPointcloudData $seeds \
--scalar microstructure,weight \
--data microstructure,weight \
--inplace \
--vtk ${seeds%.*}.vtp \