some small iprovements: reading in geometry in mesh, readability of math, output formatting in fesolving and debug, hdf5 in constitutive

This commit is contained in:
Martin Diehl 2014-05-15 08:52:16 +00:00
parent 9ffcaa03e8
commit 715e7fd918
6 changed files with 22 additions and 20 deletions

View File

@ -87,7 +87,7 @@ subroutine FE_init
write(6,'(/,a)') ' <<<+- FEsolving init -+>>>' write(6,'(/,a)') ' <<<+- FEsolving init -+>>>'
write(6,'(a)') ' $Id$' write(6,'(a)') ' $Id$'
write(6,'(a16,a)') ' Current time : ',IO_timeStamp() write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
#include "compilation_info.f90" #include "compilation_info.f90"
modelName = getSolverJobName() modelName = getSolverJobName()

View File

@ -327,7 +327,7 @@ subroutine constitutive_j2_init(fileUnit)
if (any(numerics_integrator == 4_pInt)) & if (any(numerics_integrator == 4_pInt)) &
allocate(plasticState(phase)%RK4dotState (sizeDotState,NofMyPhase),source=0.0_pReal) allocate(plasticState(phase)%RK4dotState (sizeDotState,NofMyPhase),source=0.0_pReal)
if (any(numerics_integrator == 5_pInt)) & if (any(numerics_integrator == 5_pInt)) &
allocate(plasticState(phase)%RKCK45dotState (sizeDotState,NofMyPhase),source=0.0_pReal) allocate(plasticState(phase)%RKCK45dotState (6,sizeDotState,NofMyPhase),source=0.0_pReal)
#endif #endif
endif endif
enddo initializeInstances enddo initializeInstances

View File

@ -131,7 +131,7 @@ subroutine debug_init
write(6,'(/,a)') ' <<<+- debug init -+>>>' write(6,'(/,a)') ' <<<+- debug init -+>>>'
write(6,'(a)') ' $Id$' write(6,'(a)') ' $Id$'
write(6,'(a16,a)') ' Current time : ',IO_timeStamp() write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
#include "compilation_info.f90" #include "compilation_info.f90"
if (allocated(debug_StressLoopDistribution)) & if (allocated(debug_StressLoopDistribution)) &

View File

@ -1308,9 +1308,9 @@ pure function math_qRot(Q,v)
enddo enddo
enddo enddo
math_qRot(1) = -v(1)*(T(3,3)+T(4,4)) + v(2)*(T(3,2)-T(4,1)) + v(3)*(T(4,2)+T(3,1)) math_qRot = [-v(1)*(T(3,3)+T(4,4)) + v(2)*(T(3,2)-T(4,1)) + v(3)*(T(4,2)+T(3,1)), &
math_qRot(2) = v(1)*(T(3,2)+T(4,1)) - v(2)*(T(2,2)+T(4,4)) + v(3)*(T(4,3)-T(2,1)) v(1)*(T(3,2)+T(4,1)) - v(2)*(T(2,2)+T(4,4)) + v(3)*(T(4,3)-T(2,1)), &
math_qRot(3) = v(1)*(T(4,2)-T(3,1)) + v(2)*(T(4,3)+T(2,1)) - v(3)*(T(2,2)+T(3,3)) v(1)*(T(4,2)-T(3,1)) + v(2)*(T(4,3)+T(2,1)) - v(3)*(T(2,2)+T(3,3))]
math_qRot = 2.0_pReal * math_qRot + v math_qRot = 2.0_pReal * math_qRot + v
@ -1337,8 +1337,8 @@ pure function math_RtoEuler(R)
! calculate PHI ! calculate PHI
myVal=R(3,3)/sqhkl myVal=R(3,3)/sqhkl
if(myVal > 1.0_pReal) myVal = 1.0_pReal myVal = min(myVal, 1.0_pReal)
if(myVal < -1.0_pReal) myVal = -1.0_pReal myVal = max(myVal,-1.0_pReal)
math_RtoEuler(2) = acos(myVal) math_RtoEuler(2) = acos(myVal)
@ -1347,23 +1347,23 @@ pure function math_RtoEuler(R)
math_RtoEuler(3) = 0.0_pReal math_RtoEuler(3) = 0.0_pReal
! calculate phi1 ! calculate phi1
myVal=R(1,1)/squvw myVal=R(1,1)/squvw
if(myVal > 1.0_pReal) myVal = 1.0_pReal myVal = min(myVal, 1.0_pReal)
if(myVal < -1.0_pReal) myVal = -1.0_pReal myVal = max(myVal,-1.0_pReal)
math_RtoEuler(1) = acos(myVal) math_RtoEuler(1) = acos(myVal)
if(R(2,1) > 0.0_pReal) math_RtoEuler(1) = 2.0_pReal*pi-math_RtoEuler(1) if(R(2,1) > 0.0_pReal) math_RtoEuler(1) = 2.0_pReal*pi-math_RtoEuler(1)
else else
! calculate phi2 ! calculate phi2
myVal=R(2,3)/sqhk myVal=R(2,3)/sqhk
if(myVal > 1.0_pReal) myVal = 1.0_pReal myVal = min(myVal, 1.0_pReal)
if(myVal < -1.0_pReal) myVal = -1.0_pReal myVal = max(myVal,-1.0_pReal)
math_RtoEuler(3) = acos(myVal) math_RtoEuler(3) = acos(myVal)
if(R(1,3) < 0.0) math_RtoEuler(3) = 2.0_pReal*pi-math_RtoEuler(3) if(R(1,3) < 0.0) math_RtoEuler(3) = 2.0_pReal*pi-math_RtoEuler(3)
! calculate phi1 ! calculate phi1
myVal=-R(3,2)/sin(math_RtoEuler(2)) myVal=-R(3,2)/sin(math_RtoEuler(2))
if(myVal > 1.0_pReal) myVal = 1.0_pReal myVal = min(myVal, 1.0_pReal)
if(myVal < -1.0_pReal) myVal = -1.0_pReal myVal = max(myVal,-1.0_pReal)
math_RtoEuler(1) = acos(myVal) math_RtoEuler(1) = acos(myVal)
if(R(3,1) < 0.0) math_RtoEuler(1) = 2.0_pReal*pi-math_RtoEuler(1) if(R(3,1) < 0.0) math_RtoEuler(1) = 2.0_pReal*pi-math_RtoEuler(1)

View File

@ -1251,7 +1251,8 @@ subroutine mesh_spectral_build_elements(fileUnit)
e, i, & e, i, &
headerLength = 0_pInt, & headerLength = 0_pInt, &
maxIntCount, & maxIntCount, &
homog homog, &
elemType
integer(pInt), dimension(:), allocatable :: & integer(pInt), dimension(:), allocatable :: &
microstructures microstructures
integer(pInt), dimension(3) :: & integer(pInt), dimension(3) :: &
@ -1295,16 +1296,17 @@ subroutine mesh_spectral_build_elements(fileUnit)
read(fileUnit,'(a65536)') line read(fileUnit,'(a65536)') line
enddo enddo
allocate (mesh_element (4_pInt+mesh_maxNnodes,mesh_NcpElems)); mesh_element = 0_pInt allocate (mesh_element (4_pInt+mesh_maxNnodes,mesh_NcpElems), source = 0_pInt)
allocate (microstructures (1_pInt+maxIntCount)); microstructures = 2_pInt allocate (microstructures (1_pInt+maxIntCount), source = 1_pInt)
elemType = FE_mapElemtype('C3D8R')
e = 0_pInt e = 0_pInt
do while (e < mesh_NcpElems .and. microstructures(1) > 0_pInt) ! fill expected number of elements, stop at end of data (or blank line!) do while (e < mesh_NcpElems .and. microstructures(1) > 0_pInt) ! fill expected number of elements, stop at end of data (or blank line!)
microstructures = IO_continuousIntValues(fileUnit,maxIntCount,dummyName,dummySet,0_pInt) ! get affected elements microstructures = IO_continuousIntValues(fileUnit,maxIntCount,dummyName,dummySet,0_pInt) ! get affected elements
do i = 1_pInt,microstructures(1_pInt) do i = 1_pInt,microstructures(1_pInt)
e = e+1_pInt ! valid element entry e = e+1_pInt ! valid element entry
mesh_element( 1,e) = e ! FE id mesh_element( 1,e) = e ! FE id
mesh_element( 2,e) = FE_mapElemtype('C3D8R') ! elem type mesh_element( 2,e) = elemType ! elem type
mesh_element( 3,e) = homog ! homogenization mesh_element( 3,e) = homog ! homogenization
mesh_element( 4,e) = microstructures(1_pInt+i) ! microstructure mesh_element( 4,e) = microstructures(1_pInt+i) ! microstructure
mesh_element( 5,e) = e + (e-1_pInt)/grid(1) + & mesh_element( 5,e) = e + (e-1_pInt)/grid(1) + &

View File

@ -71,8 +71,8 @@ module prec
previousDotState2, & previousDotState2, &
dotState_backup, & dotState_backup, &
RK4dotState, & RK4dotState, &
aTolState, & aTolState
RKCK45dotState real(pReal), pointer, dimension(:,:,:) :: RKCK45dotState
end type end type
#endif #endif