223031012
/
DAMASK_EICMD
mirror of https://github.com/achalhp/DAMASK_EICMD.git
1
0
Fork 0
Code Issues Projects Releases Wiki Activity

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

This commit is contained in:
Martin Diehl 2016-10-22 12:02:09 +02:00
parent aab3d9af31 20985f83b5
commit 4a305e72b1
19 changed files with 966 additions and 149 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
VERSION
View File

@ -1 +1 @@
v2.0.1-201-g7278605 v2.0.1-249-g54c0ba5

53
code/CPFEM.f90
View File

@ -10,7 +10,6 @@ module CPFEM
implicit none implicit none
private private
#if defined(Marc4DAMASK) || defined(Abaqus)
real(pReal), parameter, private :: & real(pReal), parameter, private :: &
CPFEM_odd_stress = 1e15_pReal, & !< return value for stress in case of ping pong dummy cycle CPFEM_odd_stress = 1e15_pReal, & !< return value for stress in case of ping pong dummy cycle
CPFEM_odd_jacobian = 1e50_pReal !< return value for jacobian in case of ping pong dummy cycle CPFEM_odd_jacobian = 1e50_pReal !< return value for jacobian in case of ping pong dummy cycle
@ -20,7 +19,6 @@ module CPFEM
CPFEM_dcsdE !< Cauchy stress tangent CPFEM_dcsdE !< Cauchy stress tangent
real(pReal), dimension (:,:,:,:), allocatable, private :: & real(pReal), dimension (:,:,:,:), allocatable, private :: &
CPFEM_dcsdE_knownGood !< known good tangent CPFEM_dcsdE_knownGood !< known good tangent
#endif
integer(pInt), public :: & integer(pInt), public :: &
cycleCounter = 0_pInt, & !< needs description cycleCounter = 0_pInt, & !< needs description
theInc = -1_pInt, & !< needs description theInc = -1_pInt, & !< needs description
@ -83,10 +81,6 @@ subroutine CPFEM_initAll(el,ip)
use IO, only: & use IO, only: &
IO_init IO_init
use DAMASK_interface use DAMASK_interface
#ifdef FEM
use FEZoo, only: &
FEZoo_init
#endif
implicit none implicit none
integer(pInt), intent(in) :: el, & !< FE el number integer(pInt), intent(in) :: el, & !< FE el number
@ -97,9 +91,6 @@ subroutine CPFEM_initAll(el,ip)
call DAMASK_interface_init ! Spectral and FEM interface to commandline call DAMASK_interface_init ! Spectral and FEM interface to commandline
call prec_init call prec_init
call IO_init call IO_init
#ifdef FEM
call FEZoo_init
#endif
call numerics_init call numerics_init
call debug_init call debug_init
call math_init call math_init
@ -138,16 +129,12 @@ subroutine CPFEM_init
debug_levelBasic, & debug_levelBasic, &
debug_levelExtensive debug_levelExtensive
use FEsolving, only: & use FEsolving, only: &
#if defined(Marc4DAMASK) || defined(Abaqus)
symmetricSolver, & symmetricSolver, &
#endif
restartRead, & restartRead, &
modelName modelName
#if defined(Marc4DAMASK) || defined(Abaqus)
use mesh, only: & use mesh, only: &
mesh_NcpElems, & mesh_NcpElems, &
mesh_maxNips mesh_maxNips
#endif
use material, only: & use material, only: &
material_phase, & material_phase, &
homogState, & homogState, &
@ -173,12 +160,10 @@ subroutine CPFEM_init
#include "compilation_info.f90" #include "compilation_info.f90"
endif mainProcess endif mainProcess
#if defined(Marc4DAMASK) || defined(Abaqus)
! initialize stress and jacobian to zero ! initialize stress and jacobian to zero
allocate(CPFEM_cs(6,mesh_maxNips,mesh_NcpElems)) ; CPFEM_cs = 0.0_pReal allocate(CPFEM_cs(6,mesh_maxNips,mesh_NcpElems)) ; CPFEM_cs = 0.0_pReal
allocate(CPFEM_dcsdE(6,6,mesh_maxNips,mesh_NcpElems)) ; CPFEM_dcsdE = 0.0_pReal allocate(CPFEM_dcsdE(6,6,mesh_maxNips,mesh_NcpElems)) ; CPFEM_dcsdE = 0.0_pReal
allocate(CPFEM_dcsdE_knownGood(6,6,mesh_maxNips,mesh_NcpElems)) ; CPFEM_dcsdE_knownGood = 0.0_pReal allocate(CPFEM_dcsdE_knownGood(6,6,mesh_maxNips,mesh_NcpElems)) ; CPFEM_dcsdE_knownGood = 0.0_pReal
#endif
! *** restore the last converged values of each essential variable from the binary file ! *** restore the last converged values of each essential variable from the binary file
if (restartRead) then if (restartRead) then
@ -243,21 +228,17 @@ subroutine CPFEM_init
enddo readHomogInstances enddo readHomogInstances
close (777) close (777)
#if defined(Marc4DAMASK) || defined(Abaqus)
call IO_read_realFile(777,'convergeddcsdE',modelName,size(CPFEM_dcsdE)) call IO_read_realFile(777,'convergeddcsdE',modelName,size(CPFEM_dcsdE))
read (777,rec=1) CPFEM_dcsdE read (777,rec=1) CPFEM_dcsdE
close (777) close (777)
#endif
restartRead = .false. restartRead = .false.
endif endif
#if defined(Marc4DAMASK) || defined(Abaqus)
if (iand(debug_level(debug_CPFEM), debug_levelBasic) /= 0) then if (iand(debug_level(debug_CPFEM), debug_levelBasic) /= 0) then
write(6,'(a32,1x,6(i8,1x))') 'CPFEM_cs: ', shape(CPFEM_cs) write(6,'(a32,1x,6(i8,1x))') 'CPFEM_cs: ', shape(CPFEM_cs)
write(6,'(a32,1x,6(i8,1x))') 'CPFEM_dcsdE: ', shape(CPFEM_dcsdE) write(6,'(a32,1x,6(i8,1x))') 'CPFEM_dcsdE: ', shape(CPFEM_dcsdE)
write(6,'(a32,1x,6(i8,1x),/)') 'CPFEM_dcsdE_knownGood: ', shape(CPFEM_dcsdE_knownGood) write(6,'(a32,1x,6(i8,1x),/)') 'CPFEM_dcsdE_knownGood: ', shape(CPFEM_dcsdE_knownGood)
write(6,'(a32,l1)') 'symmetricSolver: ', symmetricSolver write(6,'(a32,l1)') 'symmetricSolver: ', symmetricSolver
endif endif
#endif
flush(6) flush(6)
end subroutine CPFEM_init end subroutine CPFEM_init
@ -266,11 +247,7 @@ end subroutine CPFEM_init
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief perform initialization at first call, update variables and call the actual material model !> @brief perform initialization at first call, update variables and call the actual material model
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
#if defined(Marc4DAMASK) || defined(Abaqus)
subroutine CPFEM_general(mode, parallelExecution, ffn, ffn1, temperature_inp, dt, elFE, ip, cauchyStress, jacobian) subroutine CPFEM_general(mode, parallelExecution, ffn, ffn1, temperature_inp, dt, elFE, ip, cauchyStress, jacobian)
#else
subroutine CPFEM_general(mode, ffn, ffn1, dt, elFE, ip)
#endif
use numerics, only: & use numerics, only: &
defgradTolerance, & defgradTolerance, &
iJacoStiffness, & iJacoStiffness, &
@ -281,7 +258,6 @@ subroutine CPFEM_general(mode, ffn, ffn1, dt, elFE, ip)
debug_levelBasic, & debug_levelBasic, &
debug_levelExtensive, & debug_levelExtensive, &
debug_levelSelective, & debug_levelSelective, &
#if defined(Marc4DAMASK) || defined(Abaqus)
debug_stressMaxLocation, & debug_stressMaxLocation, &
debug_stressMinLocation, & debug_stressMinLocation, &
debug_jacobianMaxLocation, & debug_jacobianMaxLocation, &
@ -290,7 +266,6 @@ subroutine CPFEM_general(mode, ffn, ffn1, dt, elFE, ip)
debug_stressMin, & debug_stressMin, &
debug_jacobianMax, & debug_jacobianMax, &
debug_jacobianMin, & debug_jacobianMin, &
#endif
debug_e, & debug_e, &
debug_i debug_i
use FEsolving, only: & use FEsolving, only: &
@ -348,12 +323,10 @@ subroutine CPFEM_general(mode, ffn, ffn1, dt, elFE, ip)
use homogenization, only: & use homogenization, only: &
materialpoint_F, & materialpoint_F, &
materialpoint_F0, & materialpoint_F0, &
#if defined(Marc4DAMASK) || defined(Abaqus)
materialpoint_P, & materialpoint_P, &
materialpoint_dPdF, & materialpoint_dPdF, &
materialpoint_results, & materialpoint_results, &
materialpoint_sizeResults, & materialpoint_sizeResults, &
#endif
materialpoint_stressAndItsTangent, & materialpoint_stressAndItsTangent, &
materialpoint_postResults materialpoint_postResults
use IO, only: & use IO, only: &
@ -368,7 +341,6 @@ subroutine CPFEM_general(mode, ffn, ffn1, dt, elFE, ip)
real(pReal), dimension (3,3), intent(in) :: ffn, & !< deformation gradient for t=t0 real(pReal), dimension (3,3), intent(in) :: ffn, & !< deformation gradient for t=t0
ffn1 !< deformation gradient for t=t1 ffn1 !< deformation gradient for t=t1
integer(pInt), intent(in) :: mode !< computation mode 1: regular computation plus aging of results integer(pInt), intent(in) :: mode !< computation mode 1: regular computation plus aging of results
#if defined(Marc4DAMASK) || defined(Abaqus)
real(pReal), intent(in) :: temperature_inp !< temperature real(pReal), intent(in) :: temperature_inp !< temperature
logical, intent(in) :: parallelExecution !< flag indicating parallel computation of requested IPs logical, intent(in) :: parallelExecution !< flag indicating parallel computation of requested IPs
real(pReal), dimension(6), intent(out) :: cauchyStress !< stress vector in Mandel notation real(pReal), dimension(6), intent(out) :: cauchyStress !< stress vector in Mandel notation
@ -381,20 +353,13 @@ subroutine CPFEM_general(mode, ffn, ffn1, dt, elFE, ip)
real(pReal), dimension (3,3,3,3) :: H_sym, & real(pReal), dimension (3,3,3,3) :: H_sym, &
H, & H, &
jacobian3333 ! jacobian in Matrix notation jacobian3333 ! jacobian in Matrix notation
#else
logical, parameter :: parallelExecution = .true.
#endif
integer(pInt) elCP, & ! crystal plasticity element number integer(pInt) elCP, & ! crystal plasticity element number
i, j, k, l, m, n, ph, homog, mySource i, j, k, l, m, n, ph, homog, mySource
logical updateJaco ! flag indicating if JAcobian has to be updated logical updateJaco ! flag indicating if JAcobian has to be updated
character(len=1024) :: rankStr character(len=1024) :: rankStr
#if defined(Marc4DAMASK) || defined(Abaqus)
elCP = mesh_FEasCP('elem',elFE) elCP = mesh_FEasCP('elem',elFE)
#else
elCP = elFE
#endif
if (iand(debug_level(debug_CPFEM), debug_levelBasic) /= 0_pInt & if (iand(debug_level(debug_CPFEM), debug_levelBasic) /= 0_pInt &
.and. elCP == debug_e .and. ip == debug_i) then .and. elCP == debug_e .and. ip == debug_i) then
@ -411,13 +376,10 @@ subroutine CPFEM_general(mode, ffn, ffn1, dt, elFE, ip)
write(6,'(a,/)') '#############################################'; flush (6) write(6,'(a,/)') '#############################################'; flush (6)
endif endif
#if defined(Marc4DAMASK) || defined(Abaqus)
if (iand(mode, CPFEM_BACKUPJACOBIAN) /= 0_pInt) & if (iand(mode, CPFEM_BACKUPJACOBIAN) /= 0_pInt) &
CPFEM_dcsde_knownGood = CPFEM_dcsde CPFEM_dcsde_knownGood = CPFEM_dcsde
if (iand(mode, CPFEM_RESTOREJACOBIAN) /= 0_pInt) & if (iand(mode, CPFEM_RESTOREJACOBIAN) /= 0_pInt) &
CPFEM_dcsde = CPFEM_dcsde_knownGood CPFEM_dcsde = CPFEM_dcsde_knownGood
#endif
!*** age results and write restart data if requested !*** age results and write restart data if requested
if (iand(mode, CPFEM_AGERESULTS) /= 0_pInt) then if (iand(mode, CPFEM_AGERESULTS) /= 0_pInt) then
@ -514,11 +476,9 @@ subroutine CPFEM_general(mode, ffn, ffn1, dt, elFE, ip)
enddo writeHomogInstances enddo writeHomogInstances
close (777) close (777)
#if defined(Marc4DAMASK) || defined(Abaqus)
call IO_write_jobRealFile(777,'convergeddcsdE',size(CPFEM_dcsdE)) call IO_write_jobRealFile(777,'convergeddcsdE',size(CPFEM_dcsdE))
write (777,rec=1) CPFEM_dcsdE write (777,rec=1) CPFEM_dcsdE
close (777) close (777)
#endif
endif endif
endif ! results aging endif ! results aging
@ -529,22 +489,18 @@ subroutine CPFEM_general(mode, ffn, ffn1, dt, elFE, ip)
!* If no parallel execution is required, there is no need to collect FEM input !* If no parallel execution is required, there is no need to collect FEM input
if (.not. parallelExecution) then if (.not. parallelExecution) then
#if defined(Marc4DAMASK) || defined(Abaqus)
temperature(material_homog(ip,elCP))%p(thermalMapping(material_homog(ip,elCP))%p(ip,elCP)) = & temperature(material_homog(ip,elCP))%p(thermalMapping(material_homog(ip,elCP))%p(ip,elCP)) = &
temperature_inp temperature_inp
#endif
materialpoint_F0(1:3,1:3,ip,elCP) = ffn materialpoint_F0(1:3,1:3,ip,elCP) = ffn
materialpoint_F(1:3,1:3,ip,elCP) = ffn1 materialpoint_F(1:3,1:3,ip,elCP) = ffn1
elseif (iand(mode, CPFEM_COLLECT) /= 0_pInt) then elseif (iand(mode, CPFEM_COLLECT) /= 0_pInt) then
#if defined(Marc4DAMASK) || defined(Abaqus)
call random_number(rnd) call random_number(rnd)
if (rnd < 0.5_pReal) rnd = rnd - 1.0_pReal if (rnd < 0.5_pReal) rnd = rnd - 1.0_pReal
CPFEM_cs(1:6,ip,elCP) = rnd * CPFEM_odd_stress CPFEM_cs(1:6,ip,elCP) = rnd * CPFEM_odd_stress
CPFEM_dcsde(1:6,1:6,ip,elCP) = CPFEM_odd_jacobian * math_identity2nd(6) CPFEM_dcsde(1:6,1:6,ip,elCP) = CPFEM_odd_jacobian * math_identity2nd(6)
temperature(material_homog(ip,elCP))%p(thermalMapping(material_homog(ip,elCP))%p(ip,elCP)) = & temperature(material_homog(ip,elCP))%p(thermalMapping(material_homog(ip,elCP))%p(ip,elCP)) = &
temperature_inp temperature_inp
#endif
materialpoint_F0(1:3,1:3,ip,elCP) = ffn materialpoint_F0(1:3,1:3,ip,elCP) = ffn
materialpoint_F(1:3,1:3,ip,elCP) = ffn1 materialpoint_F(1:3,1:3,ip,elCP) = ffn1
CPFEM_calc_done = .false. CPFEM_calc_done = .false.
@ -569,12 +525,10 @@ subroutine CPFEM_general(mode, ffn, ffn1, dt, elFE, ip)
endif endif
outdatedFFN1 = .true. outdatedFFN1 = .true.
endif endif
#if defined(Marc4DAMASK) || defined(Abaqus)
call random_number(rnd) call random_number(rnd)
if (rnd < 0.5_pReal) rnd = rnd - 1.0_pReal if (rnd < 0.5_pReal) rnd = rnd - 1.0_pReal
CPFEM_cs(1:6,ip,elCP) = rnd*CPFEM_odd_stress CPFEM_cs(1:6,ip,elCP) = rnd*CPFEM_odd_stress
CPFEM_dcsde(1:6,1:6,ip,elCP) = CPFEM_odd_jacobian*math_identity2nd(6) CPFEM_dcsde(1:6,1:6,ip,elCP) = CPFEM_odd_jacobian*math_identity2nd(6)
#endif
!*** deformation gradient is not outdated !*** deformation gradient is not outdated
@ -607,7 +561,6 @@ subroutine CPFEM_general(mode, ffn, ffn1, dt, elFE, ip)
endif endif
!* map stress and stiffness (or return odd values if terminally ill) !* map stress and stiffness (or return odd values if terminally ill)
#if defined(Marc4DAMASK) || defined(Abaqus)
terminalIllness: if ( terminallyIll ) then terminalIllness: if ( terminallyIll ) then
call random_number(rnd) call random_number(rnd)
@ -648,11 +601,8 @@ subroutine CPFEM_general(mode, ffn, ffn1, dt, elFE, ip)
CPFEM_dcsde(1:6,1:6,ip,elCP) = math_Mandel3333to66(J_inverse * H_sym) CPFEM_dcsde(1:6,1:6,ip,elCP) = math_Mandel3333to66(J_inverse * H_sym)
endif terminalIllness endif terminalIllness
#endif
endif validCalculation endif validCalculation
#if defined(Marc4DAMASK) || defined(Abaqus)
!* report stress and stiffness !* report stress and stiffness
if ((iand(debug_level(debug_CPFEM), debug_levelExtensive) /= 0_pInt) & if ((iand(debug_level(debug_CPFEM), debug_levelExtensive) /= 0_pInt) &
.and. ((debug_e == elCP .and. debug_i == ip) & .and. ((debug_e == elCP .and. debug_i == ip) &
@ -663,11 +613,9 @@ subroutine CPFEM_general(mode, ffn, ffn1, dt, elFE, ip)
'<< CPFEM >> Jacobian/GPa at elFE ip ', elFE, ip, transpose(CPFEM_dcsdE(1:6,1:6,ip,elCP))*1.0e-9_pReal '<< CPFEM >> Jacobian/GPa at elFE ip ', elFE, ip, transpose(CPFEM_dcsdE(1:6,1:6,ip,elCP))*1.0e-9_pReal
flush(6) flush(6)
endif endif
#endif
endif endif
#if defined(Marc4DAMASK) || defined(Abaqus)
!*** warn if stiffness close to zero !*** warn if stiffness close to zero
if (all(abs(CPFEM_dcsdE(1:6,1:6,ip,elCP)) < 1e-10_pReal)) call IO_warning(601,elCP,ip) if (all(abs(CPFEM_dcsdE(1:6,1:6,ip,elCP)) < 1e-10_pReal)) call IO_warning(601,elCP,ip)
@ -696,7 +644,6 @@ subroutine CPFEM_general(mode, ffn, ffn1, dt, elFE, ip)
debug_jacobianMinLocation = [elCP, ip] debug_jacobianMinLocation = [elCP, ip]
debug_jacobianMin = minval(jacobian3333) debug_jacobianMin = minval(jacobian3333)
endif endif
#endif
end subroutine CPFEM_general end subroutine CPFEM_general

1
code/DAMASK_abaqus_std.f
View File

@ -116,7 +116,6 @@ subroutine UMAT(STRESS,STATEV,DDSDDE,SSE,SPD,SCD,&
CPFEM_CALCRESULTS, & CPFEM_CALCRESULTS, &
CPFEM_AGERESULTS, & CPFEM_AGERESULTS, &
CPFEM_COLLECT, & CPFEM_COLLECT, &
CPFEM_RESTOREJACOBIAN, &
CPFEM_BACKUPJACOBIAN, & CPFEM_BACKUPJACOBIAN, &
cycleCounter, & cycleCounter, &
theInc, & theInc, &

2
code/Makefile
View File

@ -65,7 +65,7 @@ endif
# settings for shared memory multicore support # settings for shared memory multicore support
ifeq "$(OPENMP)" "ON" ifeq "$(OPENMP)" "ON"
OPENMP_FLAG_ifort =-openmp -openmp-report0 -parallel OPENMP_FLAG_ifort =-qopenmp -parallel
OPENMP_FLAG_gfortran =-fopenmp OPENMP_FLAG_gfortran =-fopenmp
endif endif

6
code/crystallite.f90
View File

@ -513,8 +513,6 @@ subroutine crystallite_stressAndItsTangent(updateJaco)
subStepMinCryst, & subStepMinCryst, &
subStepSizeCryst, & subStepSizeCryst, &
stepIncreaseCryst, & stepIncreaseCryst, &
pert_Fg, &
pert_method, &
nCryst, & nCryst, &
numerics_integrator, & numerics_integrator, &
numerics_integrationMode, & numerics_integrationMode, &
@ -574,7 +572,6 @@ subroutine crystallite_stressAndItsTangent(updateJaco)
logical, intent(in) :: & logical, intent(in) :: &
updateJaco !< whether to update the Jacobian (stiffness) or not updateJaco !< whether to update the Jacobian (stiffness) or not
real(pReal) :: & real(pReal) :: &
myPert, & ! perturbation with correct sign
formerSubStep, & formerSubStep, &
subFracIntermediate subFracIntermediate
real(pReal), dimension(3,3) :: & real(pReal), dimension(3,3) :: &
@ -586,14 +583,11 @@ subroutine crystallite_stressAndItsTangent(updateJaco)
c, & !< counter in integration point component loop c, & !< counter in integration point component loop
i, & !< counter in integration point loop i, & !< counter in integration point loop
e, & !< counter in element loop e, & !< counter in element loop
k, &
l, &
n, startIP, endIP, & n, startIP, endIP, &
neighboring_e, & neighboring_e, &
neighboring_i, & neighboring_i, &
o, & o, &
p, & p, &
perturbation , & ! loop counter for forward,backward perturbation mode
myNcomponents, & myNcomponents, &
mySource mySource
! local variables used for calculating analytic Jacobian ! local variables used for calculating analytic Jacobian

2
code/numerics.f90
View File

@ -110,7 +110,7 @@ module numerics
fftw_plan_mode = 'FFTW_PATIENT' !< reads the planing-rigor flag, see manual on www.fftw.org, Default FFTW_PATIENT: use patient planner flag fftw_plan_mode = 'FFTW_PATIENT' !< reads the planing-rigor flag, see manual on www.fftw.org, Default FFTW_PATIENT: use patient planner flag
character(len=64), protected, public :: & character(len=64), protected, public :: &
spectral_solver = 'basicpetsc' , & !< spectral solution method spectral_solver = 'basicpetsc' , & !< spectral solution method
spectral_derivative = 'continuous' !< spectral filtering method spectral_derivative = 'continuous' !< spectral spatial derivative method
character(len=1024), protected, public :: & character(len=1024), protected, public :: &
petsc_defaultOptions = '-mech_snes_type ngmres & petsc_defaultOptions = '-mech_snes_type ngmres &
&-damage_snes_type ngmres & &-damage_snes_type ngmres &

68
code/source_thermal_externalheat.f90
View File

@ -1,6 +1,7 @@
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @author Pratheek Shanthraj, Max-Planck-Institut für Eisenforschung GmbH !> @author Pratheek Shanthraj, Max-Planck-Institut für Eisenforschung GmbH
!> @brief material subroutine for thermal source due to plastic dissipation !> @author Philip Eisenlohr, Michigan State University
!> @brief material subroutine for variable heat source
!> @details to be done !> @details to be done
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
module source_thermal_externalheat module source_thermal_externalheat
@ -10,24 +11,24 @@ module source_thermal_externalheat
implicit none implicit none
private private
integer(pInt), dimension(:), allocatable, public, protected :: & integer(pInt), dimension(:), allocatable, public, protected :: &
source_thermal_externalheat_sizePostResults, & !< cumulative size of post results source_thermal_externalheat_sizePostResults, & !< cumulative size of post results
source_thermal_externalheat_offset, & !< which source is my current thermal dissipation mechanism? source_thermal_externalheat_offset, & !< which source is my current thermal dissipation mechanism?
source_thermal_externalheat_instance !< instance of thermal dissipation source mechanism source_thermal_externalheat_instance !< instance of thermal dissipation source mechanism
integer(pInt), dimension(:,:), allocatable, target, public :: & integer(pInt), dimension(:,:), allocatable, target, public :: &
source_thermal_externalheat_sizePostResult !< size of each post result output source_thermal_externalheat_sizePostResult !< size of each post result output
character(len=64), dimension(:,:), allocatable, target, public :: & character(len=64), dimension(:,:), allocatable, target, public :: &
source_thermal_externalheat_output !< name of each post result output source_thermal_externalheat_output !< name of each post result output
integer(pInt), dimension(:), allocatable, target, public :: & integer(pInt), dimension(:), allocatable, target, public :: &
source_thermal_externalheat_Noutput !< number of outputs per instance of this source source_thermal_externalheat_Noutput !< number of outputs per instance of this source
integer(pInt), dimension(:), allocatable, private :: & integer(pInt), dimension(:), allocatable, private :: &
source_thermal_externalheat_nIntervals source_thermal_externalheat_nIntervals
real(pReal), dimension(:,:), allocatable, private :: & real(pReal), dimension(:,:), allocatable, private :: &
source_thermal_externalheat_time, & source_thermal_externalheat_time, &
source_thermal_externalheat_rate source_thermal_externalheat_rate
@ -136,23 +137,26 @@ subroutine source_thermal_externalheat_init(fileUnit)
if (phase > 0_pInt ) then; if (any(phase_source(:,phase) == SOURCE_thermal_externalheat_ID)) then ! do not short-circuit here (.and. with next if statemen). It's not safe in Fortran if (phase > 0_pInt ) then; if (any(phase_source(:,phase) == SOURCE_thermal_externalheat_ID)) then ! do not short-circuit here (.and. with next if statemen). It's not safe in Fortran
instance = source_thermal_externalheat_instance(phase) ! which instance of my source is present phase instance = source_thermal_externalheat_instance(phase) ! which instance of my source is present phase
chunkPos = IO_stringPos(line) chunkPos = IO_stringPos(line)
tag = IO_lc(IO_stringValue(line,chunkPos,1_pInt)) ! extract key tag = IO_lc(IO_stringValue(line,chunkPos,1_pInt)) ! extract key
select case(tag) select case(tag)
case ('externalheat_time') case ('externalheat_time','externalheat_rate')
if (chunkPos(1) <= 2_pInt) & if (chunkPos(1) <= 2_pInt) &
call IO_error(150_pInt,ext_msg=trim(tag)//' ('//SOURCE_thermal_externalheat_label//')') call IO_error(150_pInt,ext_msg=trim(tag)//' ('//SOURCE_thermal_externalheat_label//')')
if ( source_thermal_externalheat_nIntervals(instance) > 0_pInt &
.and. source_thermal_externalheat_nIntervals(instance) /= chunkPos(1) - 2_pInt) &
call IO_error(150_pInt,ext_msg=trim(tag)//' ('//SOURCE_thermal_externalheat_label//')')
source_thermal_externalheat_nIntervals(instance) = chunkPos(1) - 2_pInt source_thermal_externalheat_nIntervals(instance) = chunkPos(1) - 2_pInt
do interval = 1, source_thermal_externalheat_nIntervals(instance) + 1_pInt do interval = 1, source_thermal_externalheat_nIntervals(instance) + 1_pInt
temp_time(instance, interval) = IO_floatValue(line,chunkPos,1_pInt + interval) select case(tag)
case ('externalheat_time')
temp_time(instance, interval) = IO_floatValue(line,chunkPos,1_pInt + interval)
case ('externalheat_rate')
temp_rate(instance, interval) = IO_floatValue(line,chunkPos,1_pInt + interval)
end select
enddo enddo
case ('externalheat_rate')
do interval = 1, source_thermal_externalheat_nIntervals(instance) + 1_pInt
temp_rate(instance, interval) = IO_floatValue(line,chunkPos,1_pInt + interval)
enddo
end select end select
endif; endif endif; endif
enddo parsingFile enddo parsingFile
@ -162,13 +166,13 @@ subroutine source_thermal_externalheat_init(fileUnit)
initializeInstances: do phase = 1_pInt, material_Nphase initializeInstances: do phase = 1_pInt, material_Nphase
if (any(phase_source(:,phase) == SOURCE_thermal_externalheat_ID)) then if (any(phase_source(:,phase) == SOURCE_thermal_externalheat_ID)) then
NofMyPhase=count(material_phase==phase) NofMyPhase = count(material_phase==phase)
instance = source_thermal_externalheat_instance(phase) instance = source_thermal_externalheat_instance(phase)
sourceOffset = source_thermal_externalheat_offset(phase) sourceOffset = source_thermal_externalheat_offset(phase)
source_thermal_externalheat_time(instance,1:source_thermal_externalheat_nIntervals(instance)+1_pInt) = & source_thermal_externalheat_time(instance,1:source_thermal_externalheat_nIntervals(instance)+1_pInt) = &
temp_time(instance,1:source_thermal_externalheat_nIntervals(instance)+1_pInt) temp_time(instance,1:source_thermal_externalheat_nIntervals(instance)+1_pInt)
source_thermal_externalheat_rate(instance,1:source_thermal_externalheat_nIntervals(instance)+1_pInt) = & source_thermal_externalheat_rate(instance,1:source_thermal_externalheat_nIntervals(instance)+1_pInt) = &
temp_rate(instance,1:source_thermal_externalheat_nIntervals(instance)+1_pInt) temp_rate(instance,1:source_thermal_externalheat_nIntervals(instance)+1_pInt)
sizeDotState = 1_pInt sizeDotState = 1_pInt
sizeDeltaState = 0_pInt sizeDeltaState = 0_pInt
@ -200,7 +204,8 @@ subroutine source_thermal_externalheat_init(fileUnit)
end subroutine source_thermal_externalheat_init end subroutine source_thermal_externalheat_init
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief calculates derived quantities from state !> @brief rate of change of state
!> @details state only contains current time to linearly interpolate given heat powers
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine source_thermal_externalheat_dotState(ipc, ip, el) subroutine source_thermal_externalheat_dotState(ipc, ip, el)
use material, only: & use material, only: &
@ -221,12 +226,12 @@ subroutine source_thermal_externalheat_dotState(ipc, ip, el)
constituent = phasememberAt(ipc,ip,el) constituent = phasememberAt(ipc,ip,el)
sourceOffset = source_thermal_externalheat_offset(phase) sourceOffset = source_thermal_externalheat_offset(phase)
sourceState(phase)%p(sourceOffset)%dotState(1,constituent) = 1.0_pReal sourceState(phase)%p(sourceOffset)%dotState(1,constituent) = 1.0_pReal ! state is current time
end subroutine source_thermal_externalheat_dotState end subroutine source_thermal_externalheat_dotState
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief returns local vacancy generation rate !> @brief returns local heat generation rate
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine source_thermal_externalheat_getRateAndItsTangent(TDot, dTDot_dT, ipc, ip, el) subroutine source_thermal_externalheat_getRateAndItsTangent(TDot, dTDot_dT, ipc, ip, el)
use material, only: & use material, only: &
@ -244,21 +249,24 @@ subroutine source_thermal_externalheat_getRateAndItsTangent(TDot, dTDot_dT, ipc,
integer(pInt) :: & integer(pInt) :: &
instance, phase, constituent, sourceOffset, interval instance, phase, constituent, sourceOffset, interval
real(pReal) :: & real(pReal) :: &
norm_time frac_time
phase = phaseAt(ipc,ip,el) phase = phaseAt(ipc,ip,el)
constituent = phasememberAt(ipc,ip,el) constituent = phasememberAt(ipc,ip,el)
instance = source_thermal_externalheat_instance(phase) instance = source_thermal_externalheat_instance(phase)
sourceOffset = source_thermal_externalheat_offset(phase) sourceOffset = source_thermal_externalheat_offset(phase)
do interval = 1, source_thermal_externalheat_nIntervals(instance) do interval = 1, source_thermal_externalheat_nIntervals(instance) ! scan through all rate segments
norm_time = (sourceState(phase)%p(sourceOffset)%state(1,constituent) - & frac_time = (sourceState(phase)%p(sourceOffset)%state(1,constituent) - &
source_thermal_externalheat_time(instance,interval)) / & source_thermal_externalheat_time(instance,interval)) / &
(source_thermal_externalheat_time(instance,interval+1) - & (source_thermal_externalheat_time(instance,interval+1) - &
source_thermal_externalheat_time(instance,interval)) source_thermal_externalheat_time(instance,interval)) ! fractional time within segment
if (norm_time >= 0.0_pReal .and. norm_time < 1.0_pReal) & if ( (frac_time < 0.0_pReal .and. interval == 1) &
TDot = source_thermal_externalheat_rate(instance,interval ) * (1.0_pReal - norm_time) + & .or. (frac_time >= 1.0_pReal .and. interval == source_thermal_externalheat_nIntervals(instance)) &
source_thermal_externalheat_rate(instance,interval+1) * norm_time .or. (frac_time >= 0.0_pReal .and. frac_time < 1.0_pReal) ) &
TDot = source_thermal_externalheat_rate(instance,interval ) * (1.0_pReal - frac_time) + &
source_thermal_externalheat_rate(instance,interval+1) * frac_time ! interpolate heat rate between segment boundaries...
! ...or extrapolate if outside of bounds
enddo enddo
dTDot_dT = 0.0 dTDot_dT = 0.0

2
code/spectral_interface.f90
View File

@ -152,7 +152,7 @@ subroutine DAMASK_interface_init()
write(6,'(a)') ' Make sure the file "material.config" exists in the working' write(6,'(a)') ' Make sure the file "material.config" exists in the working'
write(6,'(a)') ' directory.' write(6,'(a)') ' directory.'
write(6,'(a)') ' For further configuration place "numerics.config"' write(6,'(a)') ' For further configuration place "numerics.config"'
write(6,'(a)')' and "numerics.config" in that directory.' write(6,'(a)')' and "debug.config" in that directory.'
write(6,'(/,a)')' --restart XX' write(6,'(/,a)')' --restart XX'
write(6,'(a)') ' Reads in total increment No. XX-1 and continues to' write(6,'(a)') ' Reads in total increment No. XX-1 and continues to'
write(6,'(a)') ' calculate total increment No. XX.' write(6,'(a)') ' calculate total increment No. XX.'

4
lib/damask/asciitable.py
View File

@ -336,6 +336,7 @@ class ASCIItable():
try: try:
idx.append(int(label)-1) # column given as integer number? idx.append(int(label)-1) # column given as integer number?
except ValueError: except ValueError:
label = label[1:-1] if label[0] == label[-1] and label[0] in ('"',"'") else label # remove outermost quotations
try: try:
idx.append(self.tags.index(label)) # locate string in label list idx.append(self.tags.index(label)) # locate string in label list
except ValueError: except ValueError:
@ -348,6 +349,7 @@ class ASCIItable():
idx = int(labels)-1 # offset for python array indexing idx = int(labels)-1 # offset for python array indexing
except ValueError: except ValueError:
try: try:
labels = labels[1:-1] if labels[0] == labels[-1] and labels[0] in ('"',"'") else labels # remove outermost quotations
idx = self.tags.index(labels) idx = self.tags.index(labels)
except ValueError: except ValueError:
try: try:
@ -380,6 +382,7 @@ class ASCIItable():
while idx+myDim < len(self.tags) and self.tags[idx+myDim].startswith("%i_"%(myDim+1)): while idx+myDim < len(self.tags) and self.tags[idx+myDim].startswith("%i_"%(myDim+1)):
myDim += 1 # add while found myDim += 1 # add while found
except ValueError: # column has string label except ValueError: # column has string label
label = label[1:-1] if label[0] == label[-1] and label[0] in ('"',"'") else label # remove outermost quotations
if label in self.tags: # can be directly found? if label in self.tags: # can be directly found?
myDim = 1 # scalar by definition myDim = 1 # scalar by definition
elif '1_'+label in self.tags: # look for first entry of possible multidim object elif '1_'+label in self.tags: # look for first entry of possible multidim object
@ -399,6 +402,7 @@ class ASCIItable():
while idx+dim < len(self.tags) and self.tags[idx+dim].startswith("%i_"%(dim+1)): while idx+dim < len(self.tags) and self.tags[idx+dim].startswith("%i_"%(dim+1)):
dim += 1 # add as long as found dim += 1 # add as long as found
except ValueError: # column has string label except ValueError: # column has string label
labels = labels[1:-1] if labels[0] == labels[-1] and labels[0] in ('"',"'") else labels # remove outermost quotations
if labels in self.tags: # can be directly found? if labels in self.tags: # can be directly found?
dim = 1 # scalar by definition dim = 1 # scalar by definition
elif '1_'+labels in self.tags: # look for first entry of possible multidim object elif '1_'+labels in self.tags: # look for first entry of possible multidim object

17
lib/damask/h5table.py
View File

@ -18,9 +18,9 @@ import xml.etree.cElementTree as ET
# on Python 2&3 # # on Python 2&3 #
# ---------------------------------------------------------------- # # ---------------------------------------------------------------- #
try: try:
test=isinstance('test', unicode) test = isinstance('test', unicode)
except(NameError): except(NameError):
unicode=str unicode = str
def lables_to_path(label, dsXMLPath=None): def lables_to_path(label, dsXMLPath=None):
@ -63,7 +63,7 @@ class H5Table(object):
------ ------
del_entry() -- Force delete attributes/group/datasets (Dangerous) del_entry() -- Force delete attributes/group/datasets (Dangerous)
get_attr() -- Return attributes if possible get_attr() -- Return attributes if possible
add_attr() -- Add NEW attributes to dataset/group (please delete old first!) add_attr() -- Add NEW attributes to dataset/group (no force overwrite)
get_data() -- Retrieve data in numpy.ndarray get_data() -- Retrieve data in numpy.ndarray
add_data() -- Add dataset to H5 file add_data() -- Add dataset to H5 file
get_cmdlog() -- Return the command used to generate the data if possible. get_cmdlog() -- Return the command used to generate the data if possible.
@ -120,6 +120,17 @@ class H5Table(object):
dataType, h5f_path = lables_to_path(feature_name, dataType, h5f_path = lables_to_path(feature_name,
dsXMLPath=self.dsXMLFile) dsXMLPath=self.dsXMLFile)
with h5py.File(self.h5f_path, 'a') as h5f: with h5py.File(self.h5f_path, 'a') as h5f:
# NOTE:
# --> If dataset exists, delete the old one so as to write
# a new one. For brand new dataset. For brand new one,
# record its state as fresh in the cmd log.
try:
del h5f[h5f_path]
print "***deleting old {} from {}".format(feature_name,
self.h5f_path)
except:
# if no cmd log, None will used
cmd_log = str(cmd_log) + " [FRESH]"
h5f.create_dataset(h5f_path, data=dataset) h5f.create_dataset(h5f_path, data=dataset)
# store the cmd in log is possible # store the cmd in log is possible
if cmd_log is not None: if cmd_log is not None:

103
processing/post/ascii2hdf5.py
View File

@ -28,49 +28,48 @@ from optparse import OptionParser
scriptName = os.path.splitext(os.path.basename(__file__))[0] scriptName = os.path.splitext(os.path.basename(__file__))[0]
scriptID = ' '.join([scriptName,damask.version]) scriptID = ' '.join([scriptName, damask.version])
# ----- helper function ----- # # ----- helper function ----- #
def get_rectMshVectors(xyz_array, posNum): def get_rectMshVectors(xyz_array, posNum):
"""take in a xyz array from rectangular mesh and figure out Vx, Vy, Vz""" """Get Vx, Vy, Vz for rectLinear grid"""
# need some improvement, and only works for rectangular grid # need some improvement, and only works for rectangular grid
v = sorted(list(set(xyz_array[:, posNum]))) v = sorted(list(set(xyz_array[:, posNum])))
v_interval = (v[2]+v[1])/2.0 - (v[1]+v[0])/2.0 v_interval = (v[2]+v[1])/2.0 - (v[1]+v[0])/2.0
v_start = (v[1]+v[0])/2.0 - v_interval v_start = (v[1]+v[0])/2.0 - v_interval
v_end = (v[-1]+v[-2])/2.0 + v_interval v_end = (v[-1]+v[-2])/2.0 + v_interval
V = np.linspace(v_start, v_end, len(v)+1) V = np.linspace(v_start, v_end, len(v)+1)
return V return V
# ----- MAIN ---- # # ----- MAIN ---- #
desp_msg = "Convert DAMASK ascii table to HDF5 file" desp_msg = "Convert DAMASK ascii table to HDF5 file"
parser = OptionParser(option_class=damask.extendableOption, parser = OptionParser(option_class=damask.extendableOption,
usage='%prog options [file[s]]', usage='%prog options [file[s]]',
description = desp_msg, description=desp_msg,
version = scriptID) version=scriptID)
parser.add_option('-D', '--DefinitionFile', parser.add_option('-D', '--DefinitionFile',
dest = 'storage definition file', dest='storage definition file',
type = 'string', type='string',
metavar = 'string', metavar='string',
help = 'definition file for H5 data storage') help='definition file for H5 data storage')
parser.add_option('-p', parser.add_option('-p', '--pos', '--position',
'--pos', '--position', dest='pos',
dest = 'pos', type='string', metavar='string',
type = 'string', metavar = 'string', help='label of coordinates [%default]')
help = 'label of coordinates [%default]')
parser.set_defaults(DefinitionFile='default', parser.set_defaults(DefinitionFile='default',
pos='pos') pos='pos')
(options,filenames) = parser.parse_args() (options, filenames) = parser.parse_args()
filename = filenames[0] filename = filenames[0]
if options.DefinitionFile == 'default': if options.DefinitionFile == 'default':
defFile = None defFile = None
else: else:
defFile = options.DefinitionFile defFile = options.DefinitionFile
# ----- read in data using DAMASK ASCII table class ----- # # ----- read in data using DAMASK ASCII table class ----- #
asciiTable = damask.ASCIItable(name=filename, buffered=False) asciiTable = damask.ASCIItable(name=filename, buffered=False)
@ -92,44 +91,52 @@ Vz = get_rectMshVectors(xyz_array, 2)
# use the dimension of the rectangular grid to reshape all other data # use the dimension of the rectangular grid to reshape all other data
mshGridDim = [len(Vx)-1, len(Vy)-1, len(Vz)-1] mshGridDim = [len(Vx)-1, len(Vy)-1, len(Vz)-1]
# ----- compose cmd log ----- #
cmd_log = " ".join([scriptID, filename])
# ----- create a new HDF5 file and save the data -----# # ----- create a new HDF5 file and save the data -----#
# force remove existing HDF5 file # force remove existing HDF5 file
h5fName = filename.replace(".txt", ".h5") h5fName = filename.replace(".txt", ".h5")
try: try:
os.remove(h5fName) os.remove(h5fName)
except OSError: except OSError:
pass pass
h5f = damask.H5Table(h5fName, h5f = damask.H5Table(h5fName,
new_file=True, new_file=True,
dsXMLFile=defFile) dsXMLFile=defFile)
# adding increment number as root level attributes # adding increment number as root level attributes
h5f.add_attr('inc', incNum) h5f.add_attr('inc', incNum)
# add the mesh grid data now # add the mesh grid data now
h5f.add_data("Vx", Vx) h5f.add_data("Vx", Vx, cmd_log=cmd_log)
h5f.add_data("Vy", Vy) h5f.add_data("Vy", Vy, cmd_log=cmd_log)
h5f.add_data("Vz", Vz) h5f.add_data("Vz", Vz, cmd_log=cmd_log)
# add the rest of data from table # add the rest of data from table
labelsProcessed = ['inc'] labelsProcessed = ['inc']
for fi in xrange(len(labels)): for fi in xrange(len(labels)):
featureName = labels[fi] featureName = labels[fi]
# remove trouble maker "("" and ")" from label/feature name # remove trouble maker "("" and ")" from label/feature name
if "(" in featureName: featureName = featureName.replace("(", "") if "(" in featureName:
if ")" in featureName: featureName = featureName.replace(")", "") featureName = featureName.replace("(", "")
# skip increment and duplicated columns in the ASCII table if ")" in featureName:
if featureName in labelsProcessed: continue featureName = featureName.replace(")", "")
# skip increment and duplicated columns in the ASCII table
if featureName in labelsProcessed:
continue
featureIdx = labels_idx[fi] featureIdx = labels_idx[fi]
featureDim = featuresDim[fi] featureDim = featuresDim[fi]
# grab the data hook # grab the data hook
dataset = fullTable[:, featureIdx:featureIdx+featureDim] dataset = fullTable[:, featureIdx:featureIdx+featureDim]
# mapping 2D data onto a 3D rectangular mesh to get 4D data # mapping 2D data onto a 3D rectangular mesh to get 4D data
# WARNING: In paraview, the data for a recmesh is mapped as: # WARNING: In paraview, the data for a recmesh is mapped as:
# --> len(z), len(y), len(x), size(data) # --> len(z), len(y), len(x), size(data)
# dataset = dataset.reshape((mshGridDim[0], mshGridDim[1], mshGridDim[2], # dataset = dataset.reshape((mshGridDim[0],
# dataset.shape[1])) # mshGridDim[1],
# write out data # mshGridDim[2],
print "adding {}...".format(featureName) # dataset.shape[1]))
h5f.add_data(featureName, dataset) # write out data
# write down the processed label print "adding {}...".format(featureName)
labelsProcessed.append(featureName) h5f.add_data(featureName, dataset, cmd_log=cmd_log)
# write down the processed label
labelsProcessed.append(featureName)

72
processing/post/h5_addCalculation.py Executable file
View File

@ -0,0 +1,72 @@
#!/usr/bin/env python2.7
# -*- coding: UTF-8 no BOM -*-
import os
# import re
# import sys
import collections
# import math
import damask
# import numpy as np
from optparse import OptionParser
scriptName = os.path.splitext(os.path.basename(__file__))[0]
scriptID = ' '.join([scriptName, damask.version])
# ----- Helper functions ----- #
def listify(x):
return x if isinstance(x, collections.Iterable) else [x]
# --------------------------------------------------------------------
# MAIN
# --------------------------------------------------------------------
usageEx = """
usage_in_details:
Column labels are tagged by '#label#' in formulas.
Use ';' for ',' in functions. Numpy is available as 'np'.
Special variables: #_row_# -- row index
Examples:
(1) magnitude of vector -- "np.linalg.norm(#vec#)"
(2) rounded root of row number -- "round(math.sqrt(#_row_#);3)"
"""
desp = "Add or alter column(s) with derived values according to "
desp += "user-defined arithmetic operation between column(s)."
parser = OptionParser(option_class=damask.extendableOption,
usage='%prog options [file[s]]' + usageEx,
description=desp,
version=scriptID)
parser.add_option('-l', '--label',
dest='labels',
action='extend', metavar='<string LIST>',
help='(list of) new column labels')
parser.add_option('-f', '--formula',
dest='formulas',
action='extend', metavar='<string LIST>',
help='(list of) formulas corresponding to labels')
parser.add_option('-c', '--condition',
dest='condition', metavar='string',
help='condition to filter rows')
parser.set_defaults(condition=None)
(options, filenames) = parser.parse_args()
# ----- parse formulas ----- #
for i in xrange(len(options.formulas)):
options.formulas[i] = options.formulas[i].replace(';', ',')
# ----- loop over input files ----- #
for name in filenames:
try:
h5f = damask.H5Table(name, new_file=False)
except:
print "!!!Cannot process {}".format(name)
continue
damask.util.report(scriptName, name)
# Note:
# --> not immediately needed, come back later

61
processing/post/h5_addCauchy.py Executable file
View File

@ -0,0 +1,61 @@
#!/usr/bin/env python2.7
# -*- coding: UTF-8 no BOM -*-
import os
import damask
import numpy as np
from optparse import OptionParser
scriptName = os.path.splitext(os.path.basename(__file__))[0]
scriptID = ' '.join([scriptName, damask.version])
def getCauchy(f, p):
"""return Cauchy stress for given f and p"""
# [Cauchy] = (1/det(F)) * [P].[F_transpose]
f = f.reshape((3, 3))
p = p.reshape((3, 3))
return 1.0/np.linalg.det(f)*np.dot(p, f.T).reshape(9)
# --------------------------------------------------------------------
# MAIN
# --------------------------------------------------------------------
desp = "Add column(s) containing Cauchy stress based on given column(s)"
desp += "of deformation gradient and first Piola--Kirchhoff stress."
parser = OptionParser(option_class=damask.extendableOption,
usage='%prog options [file[s]]',
description=desp,
version=scriptID)
parser.add_option('-f', '--defgrad',
dest='defgrad',
type='string', metavar='string',
help='heading for deformation gradient [%default]')
parser.add_option('-p', '--stress',
dest='stress',
type='string', metavar='string',
help='heading for first Piola--Kirchhoff stress [%default]')
parser.set_defaults(defgrad='f',
stress='p')
(options, filenames) = parser.parse_args()
# ----- loop over input H5 files ----- #
for name in filenames:
try:
h5f = damask.H5Table(name, new_file=False)
except:
continue
damask.util.report(scriptName, name)
# ----- read in data ----- #
f = h5f.get_data("f")
p = h5f.get_data("p")
# ----- calculate Cauchy stress ----- #
cauchy = [getCauchy(f_i, p_i) for f_i, p_i in zip(f, p)]
# ----- write to HDF5 file ----- #
cmd_log = " ".join([scriptID, name])
h5f.add_data('Cauchy', np.array(cauchy), cmd_log=cmd_log)

145
processing/post/h5_addIPFcolor.py Executable file
View File

@ -0,0 +1,145 @@
#!/usr/bin/env python2.7
# -*- coding: UTF-8 no BOM -*-
import os
import sys
import math
import damask
import numpy as np
from optparse import OptionParser
scriptName = os.path.splitext(os.path.basename(__file__))[0]
scriptID = ' '.join([scriptName, damask.version])
# TODO
# This implementation will have to iterate through the array one
# element at a time, maybe there are some other ways to make this
# faster.
# --------------------------------------------------------------------
# MAIN
# --------------------------------------------------------------------
desp = "Add RGB color value corresponding to TSL-OIM scheme for IPF."
parser = OptionParser(option_class=damask.extendableOption,
usage='%prog options [file[s]]',
description=desp,
version=scriptID)
parser.add_option('-p', '--pole',
dest='pole',
type='float', nargs=3, metavar='float float float',
help='lab frame direction for IPF [%default]')
msg = ', '.join(damask.Symmetry.lattices[1:])
parser.add_option('-s', '--symmetry',
dest='symmetry',
type='choice', choices=damask.Symmetry.lattices[1:],
metavar='string',
help='crystal symmetry [%default] {{{}}} '.format(msg))
parser.add_option('-e', '--eulers',
dest='eulers',
type='string', metavar='string',
help='Euler angles label')
parser.add_option('-d', '--degrees',
dest='degrees',
action='store_true',
help='Euler angles are given in degrees [%default]')
parser.add_option('-m', '--matrix',
dest='matrix',
type='string', metavar='string',
help='orientation matrix label')
parser.add_option('-a',
dest='a',
type='string', metavar='string',
help='crystal frame a vector label')
parser.add_option('-b',
dest='b',
type='string', metavar='string',
help='crystal frame b vector label')
parser.add_option('-c',
dest='c',
type='string', metavar='string',
help='crystal frame c vector label')
parser.add_option('-q', '--quaternion',
dest='quaternion',
type='string', metavar='string',
help='quaternion label')
parser.set_defaults(pole=(0.0, 0.0, 1.0),
symmetry=damask.Symmetry.lattices[-1],
degrees=False)
(options, filenames) = parser.parse_args()
# safe guarding to have only one orientation representation
# use dynamic typing to group a,b,c into frame
options.frame = [options.a, options.b, options.c]
input = [options.eulers is not None,
all(options.frame),
options.matrix is not None,
options.quaternion is not None]
if np.sum(input) != 1:
parser.error('needs exactly one input format.')
# select input label that was requested (active)
label_active = np.where(input)[0][0]
(label, dim, inputtype) = [(options.eulers, 3, 'eulers'),
(options.frame, [3, 3, 3], 'frame'),
(options.matrix, 9, 'matrix'),
(options.quaternion, 4, 'quaternion')][label_active]
# rescale degrees to radians
toRadians = math.pi/180.0 if options.degrees else 1.0
# only use normalized pole
pole = np.array(options.pole)
pole /= np.linalg.norm(pole)
# ----- Loop over input files ----- #
for name in filenames:
try:
h5f = damask.H5Table(name, new_file=False)
except:
continue
damask.util.report(scriptName, name)
# extract data from HDF5 file
if inputtype == 'eulers':
orieData = h5f.get_data(label)
elif inputtype == 'matrix':
orieData = h5f.get_data(label)
orieData = orieData.reshape(orieData.shape[0], 3, 3)
elif inputtype == 'frame':
vctr_a = h5f.get_data(label[0])
vctr_b = h5f.get_data(label[1])
vctr_c = h5f.get_data(label[2])
frame = np.column_stack((vctr_a, vctr_b, vctr_c))
orieData = frame.reshape(frame.shape[0], 3, 3)
elif inputtype == 'quaternion':
orieData = h5f.get_data(label)
# calculate the IPF color
rgbArrays = np.zeros((orieData.shape[0], 3))
for ci in xrange(rgbArrays.shape[0]):
if inputtype == 'eulers':
o = damask.Orientation(Eulers=np.array(orieData[ci, :])*toRadians,
symmetry=options.symmetry).reduced()
elif inputtype == 'matrix':
o = damask.Orientation(matrix=orieData[ci, :, :].transpose(),
symmetry=options.symmetry).reduced()
elif inputtype == 'frame':
o = damask.Orientation(matrix=orieData[ci, :, :],
symmetry=options.symmetry).reduced()
elif inputtype == 'quaternion':
o = damask.Orientation(quaternion=orieData[ci, :],
symmetry=options.symmetry).reduced()
rgbArrays[ci, :] = o.IPFcolor(pole)
# compose labels/headers for IPF color (RGB)
labelIPF = 'IPF_{:g}{:g}{:g}_{sym}'.format(*options.pole,
sym=options.symmetry.lower())
# compose cmd history (go with dataset)
cmd_log = scriptID + '\t' + ' '.join(sys.argv[1:])
# write data to HDF5 file
h5f.add_data(labelIPF, rgbArrays, cmd_log=cmd_log)

85
processing/post/h5_addMises.py Executable file
View File

@ -0,0 +1,85 @@
#!/usr/bin/env python2.7
# -*- coding: UTF-8 no BOM -*-
import os
import sys
import math
import damask
import numpy as np
from optparse import OptionParser
scriptName = os.path.splitext(os.path.basename(__file__))[0]
scriptID = ' '.join([scriptName, damask.version])
# ----- Helper functions ----- #
def calcMises(what, tensor):
"""calculate von Mises equivalent"""
dev = tensor - np.trace(tensor)/3.0*np.eye(3)
symdev = 0.5*(dev+dev.T)
return math.sqrt(np.sum(symdev*symdev.T) *
{
'stress': 3.0/2.0,
'strain': 2.0/3.0,
}[what.lower()])
# --------------------------------------------------------------------
# MAIN
# --------------------------------------------------------------------
desp = "Add von Mises equivalent values for symmetric part of requested"
parser = OptionParser(option_class=damask.extendableOption,
usage='%prog options [file[s]]',
description=desp,
version=scriptID)
parser.add_option('-e', '--strain',
dest='strain',
metavar='string',
help='name of dataset containing strain tensors')
parser.add_option('-s', '--stress',
dest='stress',
metavar='string',
help='name of dataset containing stress tensors')
parser.set_defaults(strain=None, stress=None)
(options, filenames) = parser.parse_args()
# ----- Loop over input files ----- #
for name in filenames:
try:
h5f = damask.H5Table(name, new_file=False)
except:
continue
damask.util.report(scriptName, name)
# TODO:
# Could use some refactoring here
if options.stress is not None:
# extract stress tensor from HDF5
tnsr = h5f.get_data(options.stress)
# calculate von Mises equivalent row by row
vmStress = np.zeros(tnsr.shape[0])
for ri in xrange(tnsr.shape[0]):
stressTnsr = tnsr[ri, :].reshape(3, 3)
vmStress[ri] = calcMises('stress', stressTnsr)
# compose label
label = "Mises{}".format(options.stress)
# prepare log info
cmd_log = scriptID + '\t' + ' '.join(sys.argv[1:])
# write data to HDF5 file
h5f.add_data(label, vmStress, cmd_log=cmd_log)
if options.strain is not None:
tnsr = h5f.get_data(options.strain)
vmStrain = np.zeros(tnsr.shape[0])
for ri in xrange(tnsr.shape[0]):
strainTnsr = tnsr[ri, :].reshape(3, 3)
vmStrain[ri] = calcMises('strain', strainTnsr)
label = "Mises{}".format(options.strain)
cmd_log = scriptID + '\t' + ' '.join(sys.argv[1:])
h5f.add_data(label, vmStrain, cmd_log=cmd_log)

156
processing/post/h5_addStrainTensors.py Executable file
View File

@ -0,0 +1,156 @@
#!/usr/bin/env python2.7
# -*- coding: UTF-8 no BOM -*-
import os
import sys
import damask
import numpy as np
from optparse import OptionParser
scriptName = os.path.splitext(os.path.basename(__file__))[0]
scriptID = ' '.join([scriptName, damask.version])
# ----- Helper functions ----- #
def operator(stretch, strain, eigenvalues):
# Albrecht Bertram: Elasticity and Plasticity of Large Deformations
# An Introduction (3rd Edition, 2012), p. 102
return {'V#ln': np.log(eigenvalues),
'U#ln': np.log(eigenvalues),
'V#Biot': (np.ones(3, 'd') - 1.0/eigenvalues),
'U#Biot': (eigenvalues - np.ones(3, 'd')),
'V#Green': (np.ones(3, 'd') - 1.0/eigenvalues/eigenvalues)*0.5,
'U#Green': (eigenvalues*eigenvalues - np.ones(3, 'd'))*0.5,
}[stretch+'#'+strain]
def calcEPS(defgrads, stretchType, strainType):
"""calculate specific type of strain tensor"""
eps = np.zeros(defgrads.shape) # initialize container
# TODO:
# this loop can use some performance boost
# (multi-threading?)
for ri in xrange(defgrads.shape[0]):
f = defgrads[ri, :].reshape(3, 3)
U, S, Vh = np.linalg.svd(f)
R = np.dot(U, Vh) # rotation of polar decomposition
if stretchType == 'U':
stretch = np.dot(np.linalg.inv(R), f) # F = RU
elif stretchType == 'V':
stretch = np.dot(f, np.linalg.inv(R)) # F = VR
# kill nasty noisy data
stretch = np.where(abs(stretch) < 1e-12, 0, stretch)
(D, V) = np.linalg.eig(stretch)
# flip principal component with negative Eigen values
neg = np.where(D < 0.0)
D[neg] *= -1.
V[:, neg] *= -1.
# check each vector for orthogonality
# --> brutal force enforcing orthogonal base
# and re-normalize
for i, eigval in enumerate(D):
if np.dot(V[:, i], V[:, (i+1) % 3]) != 0.0:
V[:, (i+1) % 3] = np.cross(V[:, (i+2) % 3], V[:, i])
V[:, (i+1) % 3] /= np.sqrt(np.dot(V[:, (i+1) % 3],
V[:, (i+1) % 3].conj()))
# calculate requested version of strain tensor
d = operator(stretchType, strainType, D)
eps[ri] = (np.dot(V, np.dot(np.diag(d), V.T)).real).reshape(9)
return eps
# --------------------------------------------------------------------
# MAIN
# --------------------------------------------------------------------
desp = "Add column(s) containing given strains based on given stretches"
parser = OptionParser(option_class=damask.extendableOption,
usage='%prog options [file[s]]',
description=desp,
version=scriptID)
msg = 'material strains based on right Cauchy-Green deformation, i.e., C and U'
parser.add_option('-u', '--right',
dest='right',
action='store_true',
help=msg)
msg = 'spatial strains based on left Cauchy--Green deformation, i.e., B and V'
parser.add_option('-v', '--left',
dest='left',
action='store_true',
help=msg)
parser.add_option('-0', '--logarithmic',
dest='logarithmic',
action='store_true',
help='calculate logarithmic strain tensor')
parser.add_option('-1', '--biot',
dest='biot',
action='store_true',
help='calculate biot strain tensor')
parser.add_option('-2', '--green',
dest='green',
action='store_true',
help='calculate green strain tensor')
# NOTE:
# It might be easier to just calculate one type of deformation gradient
# at a time.
msg = 'heading(s) of columns containing deformation tensor values'
parser.add_option('-f', '--defgrad',
dest='defgrad',
action='extend',
metavar='<string LIST>',
help=msg)
parser.set_defaults(right=False, left=False,
logarithmic=False, biot=False, green=False,
defgrad='f')
(options, filenames) = parser.parse_args()
stretches = []
strains = []
if options.right:
stretches.append('U')
if options.left:
stretches.append('V')
if options.logarithmic:
strains.append('ln')
if options.biot:
strains.append('Biot')
if options.green:
strains.append('Green')
if options.defgrad is None:
parser.error('no data column specified.')
# ----- Loop over input files ----- #
for name in filenames:
try:
h5f = damask.H5Table(name, new_file=False)
except:
continue
damask.util.report(scriptName, name)
# extract defgrads from HDF5 storage
F = h5f.get_data(options.defgrad)
# allow calculate multiple types of strain within the
# same cmd call
for stretchType in stretches:
for strainType in strains:
# calculate strain tensor for this type
eps = calcEPS(F, stretchType, strainType)
# compose labels/headers for this strain tensor
labelsStrain = strainType + stretchType
# prepare log info
cmd_log = scriptID + '\t' + ' '.join(sys.argv[1:])
# write data to HDF5 file
h5f.add_data(labelsStrain, eps, cmd_log=cmd_log)

191
processing/post/h5_vtkAddRectilinearGridData.py Executable file
View File

@ -0,0 +1,191 @@
#!/usr/bin/env python2.7
# -*- coding: UTF-8 no BOM -*-
import os
import vtk
import damask
from vtk.util import numpy_support
from optparse import OptionParser
scriptName = os.path.splitext(os.path.basename(__file__))[0]
scriptID = ' '.join([scriptName, damask.version])
# --------------------------------------------------------------------
# MAIN
# --------------------------------------------------------------------
msg = "Add scalars, vectors, and/or an RGB tuple from"
msg += "an HDF5 to existing VTK rectilinear grid (.vtr/.vtk)."
parser = OptionParser(option_class=damask.extendableOption,
usage='%prog options [file[s]]',
description=msg,
version=scriptID)
parser.add_option('--vtk',
dest='vtk',
type='string', metavar='string',
help='VTK file name')
parser.add_option('--inplace',
dest='inplace',
action='store_true',
help='modify VTK file in-place')
parser.add_option('-r', '--render',
dest='render',
action='store_true',
help='open output in VTK render window')
parser.add_option('-d', '--data',
dest='data',
action='extend', metavar='<string LIST>',
help='scalar/vector value(s) label(s)')
parser.add_option('-t', '--tensor',
dest='tensor',
action='extend', metavar='<string LIST>',
help='tensor (3x3) value label(s)')
parser.add_option('-c', '--color',
dest='color',
action='extend', metavar='<string LIST>',
help='RGB color tuple label')
parser.add_option('-m',
'--mode',
dest='mode',
metavar='string',
type='choice', choices=['cell', 'point'],
help='cell-centered or point-centered coordinates')
parser.set_defaults(data=[],
tensor=[],
color=[],
mode='cell',
inplace=False,
render=False)
(options, filenames) = parser.parse_args()
# ----- Legacy VTK format support ----- #
if os.path.splitext(options.vtk)[1] == '.vtr':
reader = vtk.vtkXMLRectilinearGridReader()
reader.SetFileName(options.vtk)
reader.Update()
rGrid = reader.GetOutput()
elif os.path.splitext(options.vtk)[1] == '.vtk':
reader = vtk.vtkGenericDataObjectReader()
reader.SetFileName(options.vtk)
reader.Update()
rGrid = reader.GetRectilinearGridOutput()
else:
parser.error('Unsupported VTK file type extension.')
Npoints = rGrid.GetNumberOfPoints()
Ncells = rGrid.GetNumberOfCells()
# ----- Summary output (Sanity Check) ----- #
msg = '{}: {} points and {} cells...'.format(options.vtk,
Npoints,
Ncells)
damask.util.croak(msg)
# ----- Read HDF5 file ----- #
# NOTE:
# --> It is possible in the future we are trying to add data
# from different increment into the same VTK file, but
# this feature is not supported for the moment.
# --> Let it fail, if the HDF5 is invalid, python interpretor
# --> should be able to catch this error.
h5f = damask.H5Table(filenames[0], new_file=False)
# ----- Process data ----- #
featureToAdd = {'data': options.data,
'tensor': options.tensor,
'color': options.color}
VTKarray = {} # store all vtkData in dict, then ship them to file
for dataType in featureToAdd.keys():
featureNames = featureToAdd[dataType]
for featureName in featureNames:
VTKtype = vtk.VTK_DOUBLE
VTKdata = h5f.get_data(featureName)
if dataType == 'color':
VTKtype = vtk.VTK_UNSIGNED_CHAR
VTKdata = (VTKdata*255).astype(int)
elif dataType == 'tensor':
# Force symmetries tensor type data
VTKdata[:, 1] = VTKdata[:, 3] = 0.5*(VTKdata[:, 1]+VTKdata[:, 3])
VTKdata[:, 2] = VTKdata[:, 6] = 0.5*(VTKdata[:, 2]+VTKdata[:, 6])
VTKdata[:, 5] = VTKdata[:, 7] = 0.5*(VTKdata[:, 5]+VTKdata[:, 7])
# use vtk build-in numpy support to add data (much faster)
# NOTE:
# --> deep copy is necessary here, otherwise memory leak could occur
VTKarray[featureName] = numpy_support.numpy_to_vtk(num_array=VTKdata,
deep=True,
array_type=VTKtype)
VTKarray[featureName].SetName(featureName)
# ----- ship data to vtkGrid ----- #
mode = options.mode
damask.util.croak('{} mode...'.format(mode))
# NOTE:
# --> For unknown reason, Paraview only recognize one
# tensor attributes per cell, thus it would be safe
# to only add one attributes as tensor.
for dataType in featureToAdd.keys():
featureNames = featureToAdd[dataType]
for featureName in featureNames:
if dataType == 'color':
if mode == 'cell':
rGrid.GetCellData().SetScalars(VTKarray[featureName])
elif mode == 'point':
rGrid.GetPointData().SetScalars(VTKarray[featureName])
elif dataType == 'tensor':
if mode == 'cell':
rGrid.GetCellData().SetTensors(VTKarray[featureName])
elif mode == 'point':
rGrid.GetPointData().SetTensors(VTKarray[featureName])
else:
if mode == 'cell':
rGrid.GetCellData().AddArray(VTKarray[featureName])
elif mode == 'point':
rGrid.GetPointData().AddArray(VTKarray[featureName])
rGrid.Modified()
if vtk.VTK_MAJOR_VERSION <= 5:
rGrid.Update()
# ----- write Grid to VTK file ----- #
writer = vtk.vtkXMLRectilinearGridWriter()
writer.SetDataModeToBinary()
writer.SetCompressorTypeToZLib()
vtkFileN = os.path.splitext(options.vtk)[0]
vtkExtsn = '.vtr' if options.inplace else '_added.vtr'
writer.SetFileName(vtkFileN+vtkExtsn)
if vtk.VTK_MAJOR_VERSION <= 5:
writer.SetInput(rGrid)
else:
writer.SetInputData(rGrid)
writer.Write()
# ----- render results from script ----- #
if options.render:
mapper = vtk.vtkDataSetMapper()
mapper.SetInputData(rGrid)
actor = vtk.vtkActor()
actor.SetMapper(mapper)
# Create the graphics structure. The renderer renders into the
# render window. The render window interactor captures mouse events
# and will perform appropriate camera or actor manipulation
# depending on the nature of the events.
ren = vtk.vtkRenderer()
renWin = vtk.vtkRenderWindow()
renWin.AddRenderer(ren)
ren.AddActor(actor)
ren.SetBackground(1, 1, 1)
renWin.SetSize(200, 200)
iren = vtk.vtkRenderWindowInteractor()
iren.SetRenderWindow(renWin)
iren.Initialize()
renWin.Render()
iren.Start()

135
processing/post/h5_vtkRectilinearGrid.py Executable file
View File

@ -0,0 +1,135 @@
#!/usr/bin/env python2.7
# -*- coding: UTF-8 no BOM -*-
# ------------------------------------------------------------------ #
# NOTE: #
# 1. It might be a good idea to separate IO and calculation. #
# 2. Some of the calculation could be useful in other situations, #
# why not build a math_util, or math_sup module that contains #
# all the useful functions. #
# ------------------------------------------------------------------ #
import os
import vtk
import numpy as np
import damask
from optparse import OptionParser
scriptName = os.path.splitext(os.path.basename(__file__))[0]
scriptID = ' '.join([scriptName, damask.version])
# ----- HELPER FUNCTION ----- #
def getMeshFromXYZ(xyzArray, mode):
"""calc Vx,Vy,Vz vectors for vtk rectangular mesh"""
# NOTE:
# --> np.unique will automatically sort the list
# --> although not exactly n(1), but since mesh dimension should
# small anyway, so this is still light weight.
dim = xyzArray.shape[1] # 2D:2, 3D:3
coords = [np.unique(xyzArray[:, i]) for i in xrange(dim)]
if mode == 'cell':
# since x, y, z might now have the same number of elements,
# we have to deal with them individually
for ri in xrange(dim):
vctr_pt = coords[ri]
vctr_cell = np.empty(len(vctr_pt)+1)
# calculate first and last end point
vctr_cell[0] = vctr_pt[0] - 0.5*abs(vctr_pt[1] - vctr_pt[0])
vctr_cell[-1] = vctr_pt[-1] + 0.5*abs(vctr_pt[-2] - vctr_pt[-1])
for cj in xrange(1, len(vctr_cell)-1):
vctr_cell[cj] = 0.5*(vctr_pt[cj-1] + vctr_pt[cj])
# update the coords
coords[ri] = vctr_cell
if dim < 3:
coords.append([0]) # expand to a 3D with 0 for z
# auxiliary description
grid = np.array(map(len, coords), 'i')
N = grid.prod() if mode == 'point' else (grid-1).prod()
return coords, grid, N
# --------------------------------------------------------------------
# MAIN
# --------------------------------------------------------------------
msg = "Create regular voxel grid from points in an ASCIItable."
parser = OptionParser(option_class=damask.extendableOption,
usage='%prog options [file[s]]',
description=msg,
version=scriptID)
parser.add_option('-m',
'--mode',
dest='mode',
metavar='string',
type='choice', choices=['cell', 'point'],
help='cell-centered or point-centered coordinates')
parser.add_option('-p',
'--pos', '--position',
dest='pos',
type='string', metavar='string',
help='label of coordinates [%default]')
parser.set_defaults(mode='cell',
pos='pos')
(options, filenames) = parser.parse_args()
# ----- loop over input files ----- #
for name in filenames:
try:
h5f = damask.H5Table(name, new_file=False)
except:
continue
damask.util.report(scriptName, name)
# ----- read xyzArray from HDF5 file ----- #
xyzArray = h5f.get_data(options.pos)
# ----- figure out size and grid ----- #
coords, grid, N = getMeshFromXYZ(xyzArray, options.mode)
# ----- process data ----- #
rGrid = vtk.vtkRectilinearGrid()
# WARNING: list expansion does not work here as these are
# just pointers for a vtk instance. Simply put,
# DON't USE
# [<VTK_CONSTRUCTOR>] * <NUM_OF_ELEMENTS>
coordArray = [vtk.vtkDoubleArray(),
vtk.vtkDoubleArray(),
vtk.vtkDoubleArray()]
rGrid.SetDimensions(*grid)
for i, points in enumerate(coords):
for point in points:
coordArray[i].InsertNextValue(point)
rGrid.SetXCoordinates(coordArray[0])
rGrid.SetYCoordinates(coordArray[1])
rGrid.SetZCoordinates(coordArray[2])
# ----- output result ----- #
dirPath = os.path.split(name)[0]
if name:
writer = vtk.vtkXMLRectilinearGridWriter()
writer.SetCompressorTypeToZLib()
writer.SetDataModeToBinary()
# getting the name is a little bit tricky
vtkFileName = os.path.splitext(os.path.split(name)[1])[0]
vtkFileName += '_{}({})'.format(options.pos, options.mode)
vtkFileName += '.' + writer.GetDefaultFileExtension()
writer.SetFileName(os.path.join(dirPath, vtkFileName))
else:
writer = vtk.vtkDataSetWriter()
writer.SetHeader('# powered by '+scriptID)
writer.WriteToOutputStringOn()
if vtk.VTK_MAJOR_VERSION <= 5:
writer.SetInput(rGrid)
else:
writer.SetInputData(rGrid)
writer.Write()

10
processing/post/vtk_addRectilinearGridData.py
View File

@ -14,10 +14,12 @@ scriptID = ' '.join([scriptName,damask.version])
# MAIN # MAIN
# -------------------------------------------------------------------- # --------------------------------------------------------------------
parser = OptionParser(option_class=damask.extendableOption, usage='%prog options [file[s]]', description = """ msg = "Add scalars, vectors, and/or an RGB tuple from"
Add scalars, vectors, and/or an RGB tuple from an ASCIItable to existing VTK rectilinear grid (.vtr/.vtk). msg += "an ASCIItable to existing VTK rectilinear grid (.vtr/.vtk)."
parser = OptionParser(option_class=damask.extendableOption,
""", version = scriptID) usage='%prog options [file[s]]',
description = msg,
version = scriptID)
parser.add_option( '--vtk', parser.add_option( '--vtk',
dest = 'vtk', dest = 'vtk',