moved all routines from postprocessingMath to math.90, renamed the module to DAMASK, changed scripts and interfaces accordingly.

polished math.f90 (mainly added _pInt/_pReal and intent(in/out))
curl_fft is still a dummy function
This commit is contained in:
Martin Diehl 2011-12-01 12:01:13 +00:00
parent 1f1046ecc5
commit ace6851389
13 changed files with 1781 additions and 2285 deletions

View File

@ -0,0 +1,63 @@
! Copyright 2011 Max-Planck-Institut für Eisenforschung GmbH
!
! This file is part of DAMASK,
! the Düsseldorf Advanced MAterial Simulation Kit.
!
! DAMASK is free software: you can redistribute it and/or modify
! it under the terms of the GNU General Public License as published by
! the Free Software Foundation, either version 3 of the License, or
! (at your option) any later version.
!
! DAMASK is distributed in the hope that it will be useful,
! but WITHOUT ANY WARRANTY; without even the implied warranty of
! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
! GNU General Public License for more details.
!
! You should have received a copy of the GNU General Public License
! along with DAMASK. If not, see <http://www.gnu.org/licenses/>.
!
!##############################################################
!* $Id: prec.f90 1033 2011-10-20 16:46:11Z MPIE\m.diehl $
!##############################################################
MODULE prec
implicit none
! *** Precision of real and integer variables for python interfacing***
integer, parameter :: pReal = selected_real_kind(8)
integer, parameter :: pInt = selected_int_kind(9) ! up to +- 1e9
real(pReal), parameter :: DAMASK_NaN = Z'7FF0000000000001'
real(pReal), parameter :: tol_math_check = 1.0e-8_pReal
END MODULE prec
MODULE debug
use prec, only: pInt
implicit none
integer(pInt), parameter :: debug_verbosity = 0_pInt
END MODULE debug
MODULE numerics
use prec, only: pInt
implicit none
real*8, parameter :: fftw_timelimit = -1.0
integer*8, parameter :: fftw_planner_flag = 32
integer(pInt), parameter :: fixedSeed = 1_pInt
END MODULE numerics
MODULE IO
CONTAINS
subroutine IO_error(error_ID,e,i,g,ext_msg)
use prec, only: pInt
implicit none
integer(pInt), intent(in) :: error_ID
integer(pInt), optional, intent(in) :: e,i,g
character(len=*), optional, intent(in) :: ext_msg
character(len=1024) msg
select case (error_ID)
case default
print*, 'Error messages not supported when interfacing to Python'
end select
end subroutine IO_error
END MODULE IO

View File

@ -46,7 +46,7 @@ program DAMASK_spectral
!******************************************************************** !********************************************************************
use DAMASK_interface use DAMASK_interface
use prec, only: pInt, pReal use prec, only: pInt, pReal, DAMASK_NaN
use IO use IO
use debug, only: spectral_debug_verbosity use debug, only: spectral_debug_verbosity
use math use math
@ -124,7 +124,6 @@ program DAMASK_spectral
real(pReal), dimension(:,:,:,:), allocatable :: xi ! wave vector field real(pReal), dimension(:,:,:,:), allocatable :: xi ! wave vector field
integer(pInt), dimension(3) :: k_s integer(pInt), dimension(3) :: k_s
integer*8, dimension(3) :: fftw_plan ! plans for fftw (forward and backward) integer*8, dimension(3) :: fftw_plan ! plans for fftw (forward and backward)
integer*8 :: fftw_flag ! planner flag for fftw
! loop variables, convergence etc. ! loop variables, convergence etc.
real(pReal) :: time = 0.0_pReal, time0 = 0.0_pReal, timeinc ! elapsed time, begin of interval, time interval real(pReal) :: time = 0.0_pReal, time0 = 0.0_pReal, timeinc ! elapsed time, begin of interval, time interval
@ -424,21 +423,8 @@ program DAMASK_spectral
call dfftw_plan_with_nthreads(DAMASK_NumThreadsInt) call dfftw_plan_with_nthreads(DAMASK_NumThreadsInt)
endif endif
#endif #endif
call dfftw_set_timelimit(fftw_timelimit) ! is not working, have to fix it in FFTW source file call dfftw_set_timelimit(fftw_timelimit)
select case(IO_lc(fftw_planner_flag)) ! setting parameters for the plan creation of FFTW. Basically a translation from fftw3.f !*************************************************************
case('estimate','fftw_estimate') ! ordered from slow execution (but fast plan creation) to fast execution
fftw_flag = 64
case('measure','fftw_measure')
fftw_flag = 0
case('patient','fftw_patient')
fftw_flag= 32
case('exhaustive','fftw_exhaustive')
fftw_flag = 8
case default
call IO_warning(warning_ID=47_pInt,ext_msg=trim(IO_lc(fftw_planner_flag)))
fftw_flag = 32
end select
!*************************************************************
! Loop over loadcases defined in the loadcase file ! Loop over loadcases defined in the loadcase file
do loadcase = 1_pInt, N_Loadcases do loadcase = 1_pInt, N_Loadcases
!************************************************************* !*************************************************************
@ -509,14 +495,14 @@ program DAMASK_spectral
wgt = 1.0_pReal/real(res(1)*res(2)*res(3), pReal) wgt = 1.0_pReal/real(res(1)*res(2)*res(3), pReal)
call dfftw_plan_many_dft_r2c(fftw_plan(1),3,(/res(1),res(2),res(3)/),9,& call dfftw_plan_many_dft_r2c(fftw_plan(1),3,(/res(1),res(2),res(3)/),9,&
workfft,(/res(1) +2_pInt,res(2),res(3)/),1,(res(1) +2_pInt)*res(2)*res(3),& workfft,(/res(1) +2_pInt,res(2),res(3)/),1,(res(1) +2_pInt)*res(2)*res(3),&
workfft,(/res(1)/2_pInt+1_pInt,res(2),res(3)/),1,(res(1)/2_pInt+1_pInt)*res(2)*res(3),fftw_flag) workfft,(/res(1)/2_pInt+1_pInt,res(2),res(3)/),1,(res(1)/2_pInt+1_pInt)*res(2)*res(3),fftw_planner_flag)
call dfftw_plan_many_dft_c2r(fftw_plan(2),3,(/res(1),res(2),res(3)/),9,& call dfftw_plan_many_dft_c2r(fftw_plan(2),3,(/res(1),res(2),res(3)/),9,&
workfft,(/res(1)/2_pInt+1_pInt,res(2),res(3)/),1,(res(1)/2_pInt+1_pInt)*res(2)*res(3),& workfft,(/res(1)/2_pInt+1_pInt,res(2),res(3)/),1,(res(1)/2_pInt+1_pInt)*res(2)*res(3),&
workfft,(/res(1) +2_pInt,res(2),res(3)/),1,(res(1) +2_pInt)*res(2)*res(3),fftw_flag) workfft,(/res(1) +2_pInt,res(2),res(3)/),1,(res(1) +2_pInt)*res(2)*res(3),fftw_planner_flag)
if (debugDivergence) & if (debugDivergence) &
call dfftw_plan_many_dft_c2r(fftw_plan(3),3,(/res(1),res(2),res(3)/),3,& call dfftw_plan_many_dft_c2r(fftw_plan(3),3,(/res(1),res(2),res(3)/),3,&
divergence,(/res(1)/2_pInt+1_pInt,res(2),res(3)/),1,(res(1)/2_pInt+1_pInt)*res(2)*res(3),& divergence,(/res(1)/2_pInt+1_pInt,res(2),res(3)/),1,(res(1)/2_pInt+1_pInt)*res(2)*res(3),&
divergence,(/res(1) +2_pInt,res(2),res(3)/),1,(res(1) +2_pInt)*res(2)*res(3),fftw_flag) divergence,(/res(1) +2_pInt,res(2),res(3)/),1,(res(1) +2_pInt)*res(2)*res(3),fftw_planner_flag)
if (debugGeneral) then if (debugGeneral) then
!$OMP CRITICAL (write2out) !$OMP CRITICAL (write2out)
write (6,*) 'FFTW initialized' write (6,*) 'FFTW initialized'

File diff suppressed because it is too large Load Diff

View File

@ -23,6 +23,7 @@ MODULE numerics
!############################################################## !##############################################################
use prec, only: pInt, pReal use prec, only: pInt, pReal
use IO, only: IO_warning
implicit none implicit none
character(len=64), parameter :: numerics_configFile = 'numerics.config' ! name of configuration file character(len=64), parameter :: numerics_configFile = 'numerics.config' ! name of configuration file
@ -69,7 +70,8 @@ real(pReal) :: relevantStrain, & ! strain
err_stress_tolrel, & ! factor to multiply with highest stress to get err_stress_tol err_stress_tolrel, & ! factor to multiply with highest stress to get err_stress_tol
fftw_timelimit, & ! sets the timelimit of plan creation for FFTW, see manual on www.fftw.org fftw_timelimit, & ! sets the timelimit of plan creation for FFTW, see manual on www.fftw.org
rotation_tol ! tolerance of rotation specified in loadcase rotation_tol ! tolerance of rotation specified in loadcase
character(len=64) :: fftw_planner_flag ! sets the planig-rigor flag, see manual on www.fftw.org character(len=64) :: fftw_planner_string ! reads the planing-rigor flag, see manual on www.fftw.org
integer*8 :: fftw_planner_flag ! conversion of fftw_planner_string to integer, basically what is usually done in the include file of fftw
logical :: memory_efficient,& ! for fast execution (pre calculation of gamma_hat) logical :: memory_efficient,& ! for fast execution (pre calculation of gamma_hat)
divergence_correction ! correct divergence calculation in fourier space divergence_correction ! correct divergence calculation in fourier space
integer(pInt) :: itmax , & ! maximum number of iterations integer(pInt) :: itmax , & ! maximum number of iterations
@ -170,7 +172,7 @@ subroutine numerics_init()
itmax = 20_pInt ! Maximum iteration number itmax = 20_pInt ! Maximum iteration number
memory_efficient = .true. ! Precalculate Gamma-operator (81 double per point) memory_efficient = .true. ! Precalculate Gamma-operator (81 double per point)
fftw_timelimit = -1.0_pReal ! no timelimit of plan creation for FFTW fftw_timelimit = -1.0_pReal ! no timelimit of plan creation for FFTW
fftw_planner_flag ='FFTW_PATIENT' fftw_planner_string ='FFTW_PATIENT'
rotation_tol = 1.0e-12 rotation_tol = 1.0e-12
divergence_correction = .true. divergence_correction = .true.
!* Random seeding parameters !* Random seeding parameters
@ -286,8 +288,8 @@ subroutine numerics_init()
memory_efficient = IO_intValue(line,positions,2) > 0_pInt memory_efficient = IO_intValue(line,positions,2) > 0_pInt
case ('fftw_timelimit') case ('fftw_timelimit')
fftw_timelimit = IO_floatValue(line,positions,2) fftw_timelimit = IO_floatValue(line,positions,2)
case ('fftw_planner_flag') case ('fftw_planner_string')
fftw_planner_flag = IO_stringValue(line,positions,2) fftw_planner_string = IO_stringValue(line,positions,2)
case ('rotation_tol') case ('rotation_tol')
rotation_tol = IO_floatValue(line,positions,2) rotation_tol = IO_floatValue(line,positions,2)
case ('divergence_correction') case ('divergence_correction')
@ -309,6 +311,19 @@ subroutine numerics_init()
!$OMP END CRITICAL (write2out) !$OMP END CRITICAL (write2out)
endif endif
select case(IO_lc(fftw_planner_string)) ! setting parameters for the plan creation of FFTW. Basically a translation from fftw3.f
case('estimate','fftw_estimate') ! ordered from slow execution (but fast plan creation) to fast execution
fftw_planner_flag = 64
case('measure','fftw_measure')
fftw_planner_flag = 0
case('patient','fftw_patient')
fftw_planner_flag= 32
case('exhaustive','fftw_exhaustive')
fftw_planner_flag = 8
case default
call IO_warning(warning_ID=47_pInt,ext_msg=trim(IO_lc(fftw_planner_string)))
fftw_planner_flag = 32
end select
! writing parameters to output file ! writing parameters to output file
!$OMP CRITICAL (write2out) !$OMP CRITICAL (write2out)
@ -360,7 +375,8 @@ subroutine numerics_init()
else else
write(6,'(a24,x,e8.1)') ' fftw_timelimit: ',fftw_timelimit write(6,'(a24,x,e8.1)') ' fftw_timelimit: ',fftw_timelimit
endif endif
write(6,'(a24,x,a)') ' fftw_planner_flag: ',trim(fftw_planner_flag) write(6,'(a24,x,a)') ' fftw_planner_string: ',trim(fftw_planner_string)
write(6,'(a24,x,i8)') ' fftw_planner_flag: ',fftw_planner_flag
write(6,'(a24,x,e8.1)') ' rotation_tol: ',rotation_tol write(6,'(a24,x,e8.1)') ' rotation_tol: ',rotation_tol
write(6,'(a24,x,L8,/)') ' divergence_correction: ',divergence_correction write(6,'(a24,x,L8,/)') ' divergence_correction: ',divergence_correction

View File

@ -30,7 +30,7 @@ integer, parameter :: pInt = selected_int_kind(9) ! up to +- 1e9
integer, parameter :: pLongInt = 8 ! should be 64bit integer, parameter :: pLongInt = 8 ! should be 64bit
real(pReal), parameter :: tol_math_check = 1.0e-8_pReal real(pReal), parameter :: tol_math_check = 1.0e-8_pReal
real(pReal), parameter :: tol_gravityNodePos = 1.0e-100_pReal real(pReal), parameter :: tol_gravityNodePos = 1.0e-100_pReal
! NaN is precistion dependent ! NaN is precision dependent
! from http://www.hpc.unimelb.edu.au/doc/f90lrm/dfum_035.html ! from http://www.hpc.unimelb.edu.au/doc/f90lrm/dfum_035.html
! copy found in documentation/Code/Fortran ! copy found in documentation/Code/Fortran
real(pReal), parameter :: DAMASK_NaN = Z'7FF0000000000001' real(pReal), parameter :: DAMASK_NaN = Z'7FF0000000000001'

View File

@ -5,7 +5,7 @@
# As it reads in the data coming from "materialpoint_results", it can be adopted to the data # As it reads in the data coming from "materialpoint_results", it can be adopted to the data
# computed using the FEM solvers. Its capable to handle elements with one IP in a regular order # computed using the FEM solvers. Its capable to handle elements with one IP in a regular order
import os,sys,threading,re,numpy,time,string,postprocessingMath import os,sys,threading,re,numpy,time,string,DAMASK
from optparse import OptionParser, OptionGroup, Option, SUPPRESS_HELP from optparse import OptionParser, OptionGroup, Option, SUPPRESS_HELP
# ----------------------------- # -----------------------------
@ -412,18 +412,16 @@ for filename in args:
dim[2] = options.unitlength dim[2] = options.unitlength
if options.undeformed: if options.undeformed:
defgrad_av = numpy.eye(3) defgrad_av = numpy.eye(3)
else: # @Martin: why do we have to reshape this data when just averaging?? else:
defgrad_av = postprocessingMath.tensor_avg(res[0],res[1],res[2],\ defgrad_av = DAMASK.math.tensor_avg(res,numpy.reshape(values[:,column['tensor'][options.defgrad]:
numpy.reshape(values[:,column['tensor'][options.defgrad]: column['tensor'][options.defgrad]+9],
column['tensor'][options.defgrad]+9], (res[0],res[1],res[2],3,3)))
(res[0],res[1],res[2],3,3)))
# @Martin: any reason for having 3 args for res but a single vector arg for dim? centroids = DAMASK.math.deformed_fft(res,dim,defgrad_av,options.scaling,
centroids = postprocessingMath.deformed_fft(res[0],res[1],res[2],dim,\ numpy.reshape(values[:,column['tensor'][options.defgrad]:
numpy.reshape(values[:,column['tensor'][options.defgrad]: column['tensor'][options.defgrad]+9],
column['tensor'][options.defgrad]+9], (res[0],res[1],res[2],3,3)))
(res[0],res[1],res[2],3,3)),defgrad_av,options.scaling) ms = DAMASK.math.mesh_regular_grid(res,dim,defgrad_av,centroids)
ms = postprocessingMath.mesh(res[0],res[1],res[2],dim,defgrad_av,centroids)
fields = {\ fields = {\
'tensor': {},\ 'tensor': {},\

160
processing/post/DAMASK.pyf Normal file
View File

@ -0,0 +1,160 @@
! $Id: postprocessingMath.pyf 979 2011-08-25 18:18:38Z MPIE\m.diehl $
! -*- f90 -*-
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! Note: the context of this file is case sensitive.
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! This file was auto-generated with f2py (version:2_5972).
! See http://cens.ioc.ee/projects/f2py2e/
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! The auto-generated file is quite heavily corrected
! For modifying, notice the following hints:
! - if the dimension of an array depend on a array that is itself an input, use the C-Syntax: (1) becomes [0] etc.
! - be sure that the precision defined for math.f90 is integer, real*8, and complex*16
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
python module DAMASK ! in
interface ! in :DAMASK
module math ! in :math:math.f90
subroutine volume_compare(res,geomdim,defgrad,nodes,volume_mismatch) ! in :math:math.f90
! input variables
integer, dimension(3), intent(in) :: res
real*8, dimension(3), intent(in) :: geomdim
real*8, dimension(res[0], res[1], res[2], 3,3),intent(in), depend(res[0],res[1],res[2]) :: defgrad
real*8, dimension(res[0]+1,res[1]+1,res[2]+1,3), intent(in), depend(res[0],res[1],res[2]) :: nodes
! output variables
real*8, dimension(res[0], res[1], res[2]), intent(out), depend(res[0],res[1],res[2])) :: volume_mismatch
end subroutine volume_compare
subroutine shape_compare(res,geomdim,defgrad,nodes,centroids,shape_mismatch) ! in :math:math.f90
! input variables
integer, dimension(3), intent(in) :: res
real*8, dimension(3), intent(in) :: geomdim
real*8, dimension(res[0], res[1], res[2], 3,3),intent(in), depend(res[0],res[1],res[2]) :: defgrad
real*8, dimension(res[0]+1,res[1]+1,res[2]+1,3), intent(in), depend(res[0],res[1],res[2]) :: nodes
real*8, dimension(res[0], res[1], res[2], 3), intent(in), depend(res[0],res[1],res[2]) :: centroids
! output variables
real*8, dimension(res[0], res[1], res[2]), intent(out), depend(res[0],res[1],res[2])) :: shape_mismatch
end subroutine shape_compare
subroutine mesh_regular_grid(res,geomdim,defgrad_av,centroids,nodes) ! in :math:math.f90
! input variables
integer, dimension(3), intent(in) :: res
real*8, dimension(3), intent(in) :: geomdim
real*8, dimension(3,3), intent(in) :: defgrad_av
real*8, dimension(res[0], res[1], res[2], 3), intent(in), depend(res[0],res[1],res[2]) :: centroids
! output variables
real*8, dimension(res[0]+1,res[1]+1,res[2]+1,3), intent(out), depend(res[0],res[1],res[2]) :: nodes
! variables with dimension depending on input
real*8, dimension(res[0]+2,res[1]+2,res[2]+2,3), depend(res[0],res[1],res[2]) :: wrappedCentroids
end subroutine mesh_regular_grid
subroutine deformed_linear(res,geomdim,defgrad_av,defgrad,coord_avgCorner) ! in :math:math.f90
! input variables
integer, dimension(3), intent(in) :: res
real*8, dimension(3), intent(in) :: geomdim
real*8, dimension(3,3), intent(in) :: defgrad_av
real*8, dimension(res[0], res[1], res[2], 3,3),intent(in), depend(res[0],res[1],res[2]) :: defgrad
! output variables
real*8, dimension(res[0], res[1], res[2], 3), intent(out), depend(res[0],res[1],res[2]) :: coord_avgCorner
! variables with dimension depending on input
real*8, dimension(8,6,res[0],res[1],res[2],3), depend(res[0],res[1],res[2]) :: coord
real*8, dimension(8,res[0],res[1],res[2],3), depend(res[0],res[1],res[2]) :: coord_avgOrder
end subroutine deformed_linear
subroutine deformed_fft(res,geomdim,defgrad_av,scaling,defgrad,coords) ! in :math:math.f90
! input variables
integer, dimension(3), intent(in) :: res
real*8, dimension(3), intent(in) :: geomdim
real*8, dimension(3,3), intent(in) :: defgrad_av
real*8, intent(in) :: scaling
real*8, dimension(res[0], res[1], res[2], 3,3),intent(in), depend(res[0],res[1],res[2]) :: defgrad
! output variables
real*8, dimension(res[0], res[1], res[2], 3), intent(out), depend(res[0],res[1],res[2]) :: coords
! variables with dimension depending on input
complex*16, dimension(res[0]/2+1,res[1],res[2],3), depend(res[0],res[1],res[2]) :: coords_fft
complex*16, dimension(res[0],res[1],res[2],3,3), depend(res[0],res[1],res[2]) :: defgrad_fft
end subroutine deformed_fft
subroutine curl_fft(res,geomdim,vec_tens,field,curl_fft) ! in :math:math.f90
! input variables
integer, dimension(3), intent(in) :: res
real*8, dimension(3), intent(in) :: geomdim
integer, intent(in) :: vec_tens
real*8, dimension(res[0], res[1], res[2], 3,vec_tens), intent(in), depend(res[0],res[1],res[2],vec_tens) :: field
! output variables
real*8, dimension(res[0], res[1], res[2], 3,vec_tens), intent(out), depend(res[0],res[1],res[2],vec_tens) :: curl_field
! variables with dimension depending on input
complex*16, dimension(res[0], res[1],res[2],3,vec_tens), depend(res[0],res[1],res[2],vec_tens) :: field_fft
complex*16, dimension(res[0]/2+1,res[1],res[2],3,vec_tens), depend(res[0],res[1],res[2],vec_tens) :: curl_field_fft
real*8, dimension(res[0]/2+1,res[1],res[2],3), depend(res[0],res[1],res[2]) :: xi
end subroutine curl_fft
subroutine divergence_fft(res,geomdim,vec_tens,field,divergence_field) ! in :math:math.f90
! input variables
integer, dimension(3), intent(in) :: res
real*8, dimension(3), intent(in) :: geomdim
integer, intent(in) :: vec_tens
real*8, dimension(res[0], res[1], res[2], 3,vec_tens), intent(in), depend(res[0],res[1],res[2],vec_tens) :: field
! output variables
real*8, dimension(res[0], res[1], res[2], vec_tens), intent(out), depend(res[0],res[1],res[2],vec_tens) :: divergence_field
! variables with dimension depending on input
complex*16, dimension(res[0], res[1],res[2],vec_tens,3),i depend(res[0],res[1],res[2],vec_tens) :: field_fft
complex*16, dimension(res[0]/2+1,res[1],res[2],vec_tens), depend(res[0],res[1],res[2],vec_tens) :: divergence_field_fft
real*8, dimension(res[0]/2+1,res[1],res[2],3), depend(res[0],res[1],res[2],3) :: xi
end subroutine divergence_fft
subroutine divergence_fdm(res,geomdim,vec_tens,order,field,divergence_field) ! in :math:math.f90
! input variables
integer dimension(3), intent(in) :: res
real*8 dimension(3), intent(in) :: geomdim
integer intent(in) :: vec_tens
integer, intent(in) :: order
real*8 dimension(res[0], res[1], res[2], 3,vec_tens), intent(in), depend(res[0],res[1],res[2],vec_tens) :: field
! output variables
real*8 dimension(res[0], res[1], res[2], vec_tens), intent(out), depend(res[0],res[1],res[2],vec_tens) :: divergence_field
end subroutine divergence_fdm
subroutine tensor_avg(res,tensor,avg) ! in :math:math.f90
! input variables
integer dimension(3), intent(in) :: res
real*8 dimension(res[0],res[1],res[2],3,3), intent(in), depend(res[0],res[1],res[2]) :: tensor
! output variables
real*8 dimension(3,3), intent(out) :: avg
end subroutine tensor_avg
subroutine logstrain_mat(res,defgrad,logstrain_field) ! in :math:math.f90
! input variables
integer, dimension(3), intent(in) :: res
real*8, dimension(res[0],res[1],res[2],3,3), intent(in), depend(res[0],res[1],res[2]) :: defgrad
! output variables
real*8, dimension(res[0],res[1],res[2],3,3), intent(out), depend(res[0],res[1],res[2]) :: logstrain_field
end subroutine logstrain_mat
subroutine logstrain_spat(res,defgrad,logstrain_field) ! in :math:math.f90
! input variables
integer, dimension(3), intent(in) :: res
real*8, dimension(res[0],res[1],res[2],3,3), intent(in), depend(res[0],res[1],res[2]) :: defgrad
! output variables
real*8, dimension(res[0],res[1],res[2],3,3), intent(out), depend(res[0],res[1],res[2]) :: logstrain_field
end subroutine logstrain_spat
subroutine calculate_cauchy(res,defgrad,p_stress,c_stress) ! in :math:math.f90
! input variables
integer, dimension(3), intent(in) :: res
real*8, dimension(res[0],res[1],res[2],3,3), intent(in), depend(res[0],res[1],res[2]) :: defgrad
real*8, dimension(res[0],res[1],res[2],3,3), intent(in), depend(res[0],res[1],res[2]) :: p_stress
! output variables
real*8, dimension(res[0],res[1],res[2],3,3), intent(out), depend(res[0],res[1],res[2]) :: c_stress
end subroutine calculate_cauchy
subroutine math_equivStrain33_field(res,tensor,vm) ! in :math:math.f90
! input variables
integer, dimension(3), intent(in) :: res
real*8, dimension(res[0],res[1],res[2],3,3),intent(in),depend(res[0],res[1],res[2]) :: tensor
! output variables
real*8, dimension(res[0],res[1],res[2]),intent(out),depend(res[0],res[1],res[2]) :: vm
end subroutine math_equivStrain33_field
end module math
end interface
end python module DAMASK

View File

@ -1,6 +1,6 @@
#!/usr/bin/python #!/usr/bin/python
import os,re,sys,math,string,numpy,postprocessingMath import os,re,sys,math,string,numpy,DAMASK
from optparse import OptionParser, Option from optparse import OptionParser, Option
# ----------------------------- # -----------------------------
@ -161,11 +161,12 @@ for file in files:
= numpy.array(map(float,items[column[datatype][label]: = numpy.array(map(float,items[column[datatype][label]:
column[datatype][label]+datainfo[datatype]['len']]),'d').reshape(3,3) column[datatype][label]+datainfo[datatype]['len']]),'d').reshape(3,3)
idx += 1 idx += 1
defgrad_av[label] = postprocessingMath.tensor_avg(options.res[0],options.res[1],options.res[2],defgrad[label]) print options.res
centroids[label] = postprocessingMath.deformed_fft(options.res[0],options.res[1],options.res[2],options.dim,defgrad[label],defgrad_av[label],1.0) defgrad_av[label] = DAMASK.math.tensor_avg(options.res,defgrad[label])
nodes[label] = postprocessingMath.mesh(options.res[0],options.res[1],options.res[2],options.dim,defgrad_av[label],centroids[label]) centroids[label] = DAMASK.math.deformed_fft(options.res,options.dim,defgrad_av[label],1.0,defgrad[label])
if options.shape: shape_mismatch[label] = postprocessingMath.shape_compare( options.res[0],options.res[1],options.res[2],options.dim,nodes[label],centroids[label],defgrad[label]) nodes[label] = DAMASK.math.mesh_regular_grid(options.res,options.dim,defgrad_av[label],centroids[label])
if options.volume: volume_mismatch[label] = postprocessingMath.volume_compare(options.res[0],options.res[1],options.res[2],options.dim,nodes[label], defgrad[label]) if options.shape: shape_mismatch[label] = DAMASK.math.shape_compare( options.res,options.dim,defgrad[label],nodes[label],centroids[label])
if options.volume: volume_mismatch[label] = DAMASK.math.volume_compare(options.res,options.dim,defgrad[label],nodes[label])
# ------------------------------------------ read file --------------------------------------- # ------------------------------------------ read file ---------------------------------------

View File

@ -1,6 +1,6 @@
#!/usr/bin/python #!/usr/bin/python
import os,re,sys,math,string,numpy,postprocessingMath import os,re,sys,math,string,numpy,DAMASK
from optparse import OptionParser, Option from optparse import OptionParser, Option
# ----------------------------- # -----------------------------
@ -21,24 +21,19 @@ class extendableOption(Option):
else: else:
Option.take_action(self, action, dest, opt, value, values, parser) Option.take_action(self, action, dest, opt, value, values, parser)
def location(idx,res): def location(idx,res):
return ( idx % res[0], \ return ( idx % res[0], \
(idx // res[0]) % res[1], \ (idx // res[0]) % res[1], \
(idx // res[0] // res[1]) % res[2] ) (idx // res[0] // res[1]) % res[2] )
def index(location,res): def index(location,res):
return ( location[0] % res[0] + \
return ( location[0] % res[0] + \ (location[1] % res[1]) * res[0] + \
(location[1] % res[1]) * res[0] + \ (location[2] % res[2]) * res[0] * res[1] )
(location[2] % res[2]) * res[0] * res[1] )
def prefixMultiply(what,len): def prefixMultiply(what,len):
return {True: ['%i_%s'%(i+1,what) for i in range(len)], return {True: ['%i_%s'%(i+1,what) for i in range(len)],
False:[what]}[len>1] False:[what]}[len>1]
# -------------------------------------------------------------------- # --------------------------------------------------------------------
@ -205,7 +200,6 @@ for file in files:
continue continue
output += '\t'.join(items) output += '\t'.join(items)
(x,y,z) = location(idx,options.res) (x,y,z) = location(idx,options.res)
for datatype,labels in active.items(): for datatype,labels in active.items():
@ -252,9 +246,9 @@ for file in files:
reshape((options.res[0],options.res[1],options.res[2],\ reshape((options.res[0],options.res[1],options.res[2],\
datainfo[datatype]['len']//3)) datainfo[datatype]['len']//3))
if accuracy == 'fft': if accuracy == 'fft':
div_field[datatype][label][accuracy] = postprocessingMath.divergence_fft(options.res[0],options.res[1],options.res[2],datainfo[datatype]['len']//3,options.dim,values[datatype][label]) div_field[datatype][label][accuracy] = DAMASK.math.divergence_fft(options.res,options.dim,datainfo[datatype]['len']//3,values[datatype][label])
else: else:
div_field[datatype][label][accuracy] = postprocessingMath.divergence(options.res[0],options.res[1],options.res[2],datainfo[datatype]['len']//3,eval(accuracy)//2-1,options.dim,values[datatype][label]) div_field[datatype][label][accuracy] = DAMASK.math.divergence_fdm(options.res,options.dim,datainfo[datatype]['len']//3,eval(accuracy)//2-1,values[datatype][label])
idx = 0 idx = 0
for line in data: for line in data:

View File

@ -0,0 +1,26 @@
#!/bin/bash
# $Id: make_postprocessingMath 1106 2011-11-17 21:36:56Z MPIE\m.diehl $
# This script is used to compile math.f90 and make the functions defined in DAMASK_math.pyf
# avialable for python in the module DAMASK_math.so
# It uses the fortran wrapper f2py that is included in the numpy package to construct the
# module postprocessingMath.so out of the fortran code postprocessingMath.f90
# for the generation of the pyf file:
#f2py -m DAMASK -h DAMASK.pyf --overwrite-signature ../../code/math.f90 \
if [[ $# -eq 0 ]]; then
wd='.'
else
wd=$1
fi
# use env. variables here!
cd $wd
rm ../../lib/DAMASK.so
f2py -c \
DAMASK.pyf \
../../code/DAMASK2Python_helper.f90 \
../../code/math.f90 \
../../lib/libfftw3.a \
/opt/acml4.4.0/ifort64/lib/libacml.a
mv DAMASK.so ../../lib/.

View File

@ -1,27 +0,0 @@
#!/bin/bash
# $Id$
# This script is used to compile the python module used for geometry reconstruction.
# It uses the fortran wrapper f2py that is included in the numpy package to construct the
# module postprocessingMath.so out of the fortran code postprocessingMath.f90
# for the generation of the pyf file:
# f2py -m postprocessingMath -h postprocessingMath.pyf --overwrite-signature postprocessingMath.f90
# use ./configure --enable-sse2 --enable-shared for the compilation of fftw3.3
if [[ $# -eq 0 ]]; then
wd='.'
else
wd=$1
fi
cd $wd
rm ../../lib/postprocessingMath.so
f2py -c \
postprocessingMath.pyf \
postprocessingMath.f90 \
../../lib/libfftw3.a \
-L./
mv postprocessingMath.so ../../lib/.

File diff suppressed because it is too large Load Diff

View File

@ -1,188 +0,0 @@
! $Id$
! -*- f90 -*-
! Note: the context of this file is case sensitive.
python module postprocessingMath ! in
interface ! in :postprocessingMath
module math ! in :postprocessingMath:postprocessingMath.f90
real*8 parameter,optional,dimension(3,3) :: math_i3=reshape((/1.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,1.0/),(/3,3/))
real*8 parameter,optional :: pi=3.14159265359
function math_mul33x33(a,b) ! in :postprocessingMath:postprocessingMath.f90:math
real*8 dimension(3,3),intent(in) :: a
real*8 dimension(3,3),intent(in) :: b
real*8 dimension(3,3) :: math_mul33x33
end function math_mul33x33
subroutine math_invert3x3(a,inva,deta,error) ! in :postprocessingMath:postprocessingMath.f90:math
real*8 dimension(3,3),intent(in) :: a
real*8 dimension(3,3),intent(out) :: inva
real*8 intent(out) :: deta
logical intent(out) :: error
end subroutine math_invert3x3
function math_det3x3(m) ! in :postprocessingMath:postprocessingMath.f90:math
real*8 dimension(3,3),intent(in) :: m
real*8 :: math_det3x3
end function math_det3x3
subroutine math_pdecomposition(fe,u,r,error) ! in :postprocessingMath:postprocessingMath.f90:math
real*8 dimension(3,3),intent(in) :: fe
real*8 dimension(3,3),intent(out) :: u
real*8 dimension(3,3),intent(out) :: r
logical intent(out) :: error
end subroutine math_pdecomposition
function math_inv3x3(a) ! in :postprocessingMath:postprocessingMath.f90:math
real*8 dimension(3,3),intent(in) :: a
real*8 dimension(3,3) :: math_inv3x3
end function math_inv3x3
subroutine math_hi(m,hi1m,hi2m,hi3m) ! in :postprocessingMath:postprocessingMath.f90:math
real*8 dimension(3,3),intent(in) :: m
real*8 intent(out) :: hi1m
real*8 intent(out) :: hi2m
real*8 intent(out) :: hi3m
end subroutine math_hi
subroutine math_spectral1(m,ew1,ew2,ew3,eb1,eb2,eb3) ! in :postprocessingMath:postprocessingMath.f90:math
real*8 dimension(3,3),intent(in) :: m
real*8 intent(out) :: ew1
real*8 intent(out) :: ew2
real*8 intent(out) :: ew3
real*8 dimension(3,3),intent(out) :: eb1
real*8 dimension(3,3),intent(out) :: eb2
real*8 dimension(3,3),intent(out) :: eb3
end subroutine math_spectral1
function math_volTetrahedron(v1,v2,v3,v4) ! in :postprocessingMath:postprocessingMath.f90:math
real*8 dimension(3), intent(in) :: v1
real*8 dimension(3), intent(in) :: v2
real*8 dimension(3), intent(in) :: v3
real*8 dimension(3), intent(in) :: v4
end function math_volTetrahedron
function mesh_location(idx,resolution) ! in :postprocessingMath:postprocessingMath.f90:math
integer, intent(in) :: idx
integer dimension(3), intent(in) :: resolution
integer dimension(3) :: mesh_location
end function mesh_location
function mesh_index(location,resolution) ! in :postprocessingMath:postprocessingMath.f90:math
integer dimension(3), intent(in) :: resolution
integer dimension(3), intent(in) :: location
integer :: mesh_index
end function mesh_location
end module math
subroutine mesh(res_x,res_y,res_z,geomdim,defgrad_av,centroids,nodes) ! in :postprocessingMath:postprocessingMath.f90
integer intent(in) :: res_x
integer intent(in) :: res_y
integer intent(in) :: res_z
real*8 dimension(3),intent(in) :: geomdim
real*8 dimension(3,3),intent(in) :: defgrad_av
real*8 dimension(res_x,res_y,res_z,3),intent(in),depend(res_x,res_y,res_z) :: centroids
real*8 dimension(res_x + 2,res_y + 2,res_z +2 ,3),depend(res_x,res_y,res_z) :: centroids
real*8 dimension(res_x + 1,res_y + 1,res_z + 1,3),intent(out),depend(res_x,res_y,res_z) :: nodes
end subroutine mesh
subroutine deformed(res_x,res_y,res_z,geomdim,defgrad,defgrad_av,coord_avgCorner) ! in :postprocessingMath:postprocessingMath.f90
integer intent(in) :: res_x
integer intent(in) :: res_y
integer intent(in) :: res_z
real*8 dimension(3),intent(in) :: geomdim
real*8 dimension(3,3),intent(in) :: defgrad_av
real*8 dimension(res_x,res_y,res_z,3),intent(out),depend(res_x,res_y,res_z) :: coord_avgCorner
real*8 dimension(8,6,res_x,res_y,res_z,3),depend(res_x,res_y,res_z) :: coord
real*8 dimension(8,res_x,res_y,res_z,3),depend(res_x,res_y,res_z) :: coord_avgOrder
real*8 dimension(res_x,res_y,res_z,3,3),intent(in),depend(res_x,res_y,res_z) :: defgrad
end subroutine deformed
subroutine deformed_fft(res_x,res_y,res_z,geomdim,defgrad,defgrad_av,scaling,coords) ! in :postprocessingMath:postprocessingMath.f90
integer intent(in) :: res_x
integer intent(in) :: res_y
integer intent(in) :: res_z
real*8 intent(in) :: scaling
real*8 dimension(3),intent(in) :: geomdim
real*8 dimension(3,3),intent(in) :: defgrad_av
real*8 dimension(res_x,res_y,res_z,3,3),intent(in),depend(res_x,res_y,res_z) :: defgrad
real*8 dimension(res_x,res_y,res_z,3),intent(out),depend(res_x,res_y,res_z) :: coords
complex*16 dimension(res_x,res_y,res_z,3,3),depend(res_x,res_y,res_z) :: defgrad_fft
complex*16 dimension(res_x/2+1,res_y,res_z,3),depend(res_x,res_y,res_z) :: coords_fft
end subroutine deformed_fft
subroutine volume_compare(res_x,res_y,res_z,geomdim,nodes,defgrad,volume_mismatch) ! in :postprocessingMath:postprocessingMath.f90
integer intent(in) :: res_x
integer intent(in) :: res_y
integer intent(in) :: res_z
real*8 dimension(3),intent(in) :: geomdim
real*8 dimension(res_x+1,res_y+1,res_z+1,3),intent(in),depend(res_x,res_y,res_z) :: nodes
real*8 dimension(res_x,res_y,res_z,3,3),intent(in),depend(res_x,res_y,res_z) :: defgrad
real*8 dimension(res_x,res_y,res_z),intent(out),depend(res_x,res_y,res_z) :: volume_mismatch
end subroutine volume_compare
subroutine shape_compare(res_x,res_y,res_z,geomdim,nodes,centroids,defgrad,shape_mismatch) ! in :postprocessingMath:postprocessingMath.f90
integer intent(in) :: res_x
integer intent(in) :: res_y
integer intent(in) :: res_z
real*8 dimension(3),intent(in) :: geomdim
real*8 dimension(res_x,res_y,res_z,3,3),intent(in),depend(res_x,res_y,res_z) :: defgrad
real*8 dimension(res_x+1,res_y+1,res_z+1,3),intent(in),depend(res_x,res_y,res_z) :: nodes
real*8 dimension(res_x,res_y,res_z,3),intent(in),depend(res_x,res_y,res_z) :: centroids
real*8 dimension(res_x,res_y,res_z),intent(out),depend(res_x,res_y,res_z) :: shape_mismatch
end subroutine shape_compare
subroutine inverse_reconstruction(res_x,res_y,res_z,reference_configuration,current_configuration,defgrad) ! in :postprocessingMath:postprocessingMath.f90
integer intent(in) :: res_x
integer intent(in) :: res_y
integer intent(in) :: res_z
real*8 dimension(res_x+1,res_y+1,res_z+1,3),intent(in),depend(res_x,res_y,res_z) :: reference_configuration
real*8 dimension(res_x+1,res_y+1,res_z+1,3),intent(in),depend(res_x,res_y,res_z) :: current_configuration
real*8 dimension(res_x,res_y,res_z,3,3),intent(out),depend(res_x,res_y,res_z) :: defgrad
end subroutine inverse_reconstruction
subroutine tensor_avg(res_x,res_y,res_z,tensor,avg) ! in :postprocessingMath:postprocessingMath.f90
integer intent(in) :: res_x
integer intent(in) :: res_y
integer intent(in) :: res_z
real*8 dimension(res_x,res_y,res_z,3,3),intent(in),depend(res_x,res_y,res_z) :: tensor
real*8 dimension(3,3),intent(out) :: avg
end subroutine tensor_avg
subroutine logstrain_mat(res_x,res_y,res_z,defgrad,logstrain_field) ! in :postprocessingMath:postprocessingMath.f90
integer intent(in) :: res_x
integer intent(in) :: res_y
integer intent(in) :: res_z
real*8 dimension(res_x,res_y,res_z,3,3),intent(in),depend(res_x,res_y,res_z) :: defgrad
real*8 dimension(res_x,res_y,res_z,3,3),intent(out),depend(res_x,res_y,res_z) :: logstrain_field
end subroutine logstrain_mat
subroutine logstrain_spat(res_x,res_y,res_z,defgrad,logstrain_field) ! in :postprocessingMath:postprocessingMath.f90
integer intent(in) :: res_x
integer intent(in) :: res_y
integer intent(in) :: res_z
real*8 dimension(res_x,res_y,res_z,3,3),intent(in),depend(res_x,res_y,res_z) :: defgrad
real*8 dimension(res_x,res_y,res_z,3,3),intent(out),depend(res_x,res_y,res_z) :: logstrain_field
end subroutine logstrain_spat
subroutine calculate_cauchy(res_x,res_y,res_z,defgrad,p_stress,c_stress) ! in :postprocessingMath:postprocessingMath.f90
integer intent(in) :: res_x
integer intent(in) :: res_y
integer intent(in) :: res_z
real*8 dimension(res_x,res_y,res_z,3,3),intent(in),depend(res_x,res_y,res_z) :: defgrad
real*8 dimension(res_x,res_y,res_z,3,3),intent(in),depend(res_x,res_y,res_z) :: p_stress
real*8 dimension(res_x,res_y,res_z,3,3),intent(out),depend(res_x,res_y,res_z) :: c_stress
end subroutine calculate_cauchy
subroutine calculate_mises(res_x,res_y,res_z,tensor,vm) ! in :postprocessingMath:postprocessingMath.f90
integer intent(in) :: res_x
integer intent(in) :: res_y
integer intent(in) :: res_z
real*8 dimension(res_x,res_y,res_z,3,3),intent(in),depend(res_x,res_y,res_z) :: tensor
real*8 dimension(res_x,res_y,res_z,1),intent(out),depend(res_x,res_y,res_z) :: vm
end subroutine calculate_mises
subroutine divergence_fft(res_x,res_y,res_z,vec_tens,geomdim,field,divergence_field) ! in :postprocessingMath:postprocessingMath.f90
integer intent(in) :: res_x
integer intent(in) :: res_y
integer intent(in) :: res_z
integer intent(in) :: vec_tens
real*8 dimension(3),intent(in) :: geomdim
real*8 dimension(res_x,res_y,res_z,3,vec_tens),intent(in),depend(res_x,res_y,res_z,vec_tens) :: field
real*8 dimension(res_x,res_y,res_z,vec_tens), intent(out),depend(res_x,res_y,res_z,vec_tens) :: divergence_field
complex*8 dimension(res_x,res_y,res_z,3,vec_tens),depend(res_x,res_y,res_z,vec_tens) :: field_fft
complex*8 dimension(res_x/2+1,res_y,res_z,vec_tens),depend(res_x,res_y,res_z,vec_tens) :: divergence_field_fft
end subroutine divergence_fft
subroutine divergence(res_x,res_y,res_z,vec_tens,order,geomdim,field,divergence_field) ! in :postprocessingMath:postprocessingMath.f90
integer intent(in) :: res_x
integer intent(in) :: res_y
integer intent(in) :: res_z
integer intent(in) :: vec_tens
integer intent(in) :: order
real*8 dimension(3),intent(in) :: geomdim
real*8 dimension(res_x,res_y,res_z,vec_tens,3),intent(in),depend(res_x,res_y,res_z,vec_tens,3) :: field
real*8 dimension(res_x,res_y,res_z,vec_tens),intent(out),depend(res_x,res_y,res_z,vec_tens) :: divergence_field
end subroutine divergence
end interface
end python module postprocessingMath
! This file was auto-generated with f2py (version:2_5972).
! See http://cens.ioc.ee/projects/f2py2e/
! modified by m.diehl