corrected bug causing wrong element construction when using spectral solver. Also improved error messages for invalid resolution, dimension, and homogenization.

This commit is contained in:
Martin Diehl 2012-05-21 09:06:02 +00:00
parent 20f203f799
commit cec3357e8f
3 changed files with 32 additions and 34 deletions

View File

@ -1302,8 +1302,6 @@ subroutine IO_error(error_ID,e,i,g,ext_msg)
!* errors related to spectral solver !* errors related to spectral solver
case (802_pInt)
msg = 'invaldid geometry parameter'
case (808_pInt) case (808_pInt)
msg = 'precision not suitable for FFTW' msg = 'precision not suitable for FFTW'
case (809_pInt) case (809_pInt)
@ -1322,8 +1320,10 @@ subroutine IO_error(error_ID,e,i,g,ext_msg)
msg = 'incomplete loadcase' msg = 'incomplete loadcase'
case (838_pInt) case (838_pInt)
msg = 'mixed boundary conditions allow rotation' msg = 'mixed boundary conditions allow rotation'
case (841_pInt)
msg = 'missing header length info in spectral mesh'
case (842_pInt) case (842_pInt)
msg = 'missing header info in spectral mesh' msg = 'homogenization in spectral mesh'
case (843_pInt) case (843_pInt)
msg = 'resolution in spectral mesh' msg = 'resolution in spectral mesh'
case (844_pInt) case (844_pInt)

View File

@ -4,7 +4,7 @@
debug # debug.f90, possible keys: basic, extensive debug # debug.f90, possible keys: basic, extensive
math # math.f90, possible key: basic math # math.f90, possible key: basic
FEsolving # FEsolving.f90, possible key: basic FEsolving # FEsolving.f90, possible key: basic
mesh # mesh.f90, possible key: basic mesh # mesh.f90, possible key: basic, extensive
material # material.f90, possible keys: basic, extensive material # material.f90, possible keys: basic, extensive
lattice # lattice.f90, possible key: basic lattice # lattice.f90, possible key: basic
constitutive # constitutive_*.f90 possible keys: basic, extensive, selective constitutive # constitutive_*.f90 possible keys: basic, extensive, selective

View File

@ -2490,21 +2490,16 @@ subroutine mesh_spectral_build_nodes(myUnit)
allocate ( mesh_node0 (3,mesh_Nnodes) ); mesh_node0 = 0.0_pReal allocate ( mesh_node0 (3,mesh_Nnodes) ); mesh_node0 = 0.0_pReal
allocate ( mesh_node (3,mesh_Nnodes) ); mesh_node = 0.0_pReal allocate ( mesh_node (3,mesh_Nnodes) ); mesh_node = 0.0_pReal
res = mesh_spectral_getResolution(myUnit) + 1_pInt res = mesh_spectral_getResolution(myUnit)
geomdim = mesh_spectral_getDimension(myUnit) geomdim = mesh_spectral_getDimension(myUnit)
if ((res(1) < 1_pInt) .or. (res(2) < 1_pInt) .or. (res(3) < 0_pInt)) &
call IO_error(error_ID=843_pInt) ! 1_pInt is already added
if ((geomdim(1) <= 0.0_pReal) .or. (geomdim(2) <= 0.0_pReal) .or. (geomdim(3) <= 0.0_pReal)) &
call IO_error(error_ID=844_pInt)
forall (n = 0_pInt:mesh_Nnodes-1_pInt) forall (n = 0_pInt:mesh_Nnodes-1_pInt)
mesh_node0(1,n+1_pInt) = geomdim(1) * real(mod(n,res(1) ),pReal) & mesh_node0(1,n+1_pInt) = geomdim(1) * real(mod(n,(res(1)+1_pInt) ),pReal) &
/ real(res(1) -1_pInt,pReal) / real(res(1),pReal)
mesh_node0(2,n+1_pInt) = geomdim(2) * real(mod(n/res(1) ,res(2)),pReal) & mesh_node0(2,n+1_pInt) = geomdim(2) * real(mod(n/(res(1)+1_pInt),(res(2)+1_pInt)),pReal) &
/ real(res(2) -1_pInt,pReal) / real(res(2),pReal)
mesh_node0(3,n+1_pInt) = geomdim(3) * real(mod(n/res(1) /res(2),res(3)),pReal) & mesh_node0(3,n+1_pInt) = geomdim(3) * real(mod(n/(res(1)+1_pInt)/(res(2)+1_pInt),(res(3)+1_pInt)),pReal) &
/ real(res(3) -1_pInt,pReal) / real(res(3),pReal)
end forall end forall
mesh_node = mesh_node0 !why? mesh_node = mesh_node0 !why?
@ -2666,7 +2661,7 @@ subroutine mesh_spectral_build_elements(myUnit)
if (keyword(1:4) == 'head') then if (keyword(1:4) == 'head') then
headerLength = IO_intValue(line,myPos,1_pInt) + 1_pInt headerLength = IO_intValue(line,myPos,1_pInt) + 1_pInt
else else
call IO_error(error_ID=842_pInt) call IO_error(error_ID=841_pInt, ext_msg='mesh_spectral_build_elements')
endif endif
rewind(myUnit) rewind(myUnit)
@ -2699,15 +2694,15 @@ subroutine mesh_spectral_build_elements(myUnit)
mesh_element( 2,e) = FE_mapElemtype('C3D8R') ! elem type mesh_element( 2,e) = FE_mapElemtype('C3D8R') ! 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)/(res(1)-1_pInt) + & mesh_element( 5,e) = e + (e-1_pInt)/res(1) + &
((e-1_pInt)/((res(1)-1_pInt)*(res(2)-1_pInt)))*res(1) ! base node ((e-1_pInt)/(res(1)*res(2)))*(res(1)+1_pInt) ! base node
mesh_element( 6,e) = mesh_element(5,e) + 1_pInt mesh_element( 6,e) = mesh_element(5,e) + 1_pInt
mesh_element( 7,e) = mesh_element(5,e) + res(1) + 1_pInt mesh_element( 7,e) = mesh_element(5,e) + res(1) + 2_pInt
mesh_element( 8,e) = mesh_element(5,e) + res(1) mesh_element( 8,e) = mesh_element(5,e) + res(1) + 1_pInt
mesh_element( 9,e) = mesh_element(5,e) + res(1) * res(2) ! second floor base node mesh_element( 9,e) = mesh_element(5,e) +(res(1) + 1_pInt) * (res(2) + 1_pInt) ! second floor base node
mesh_element(10,e) = mesh_element(9,e) + 1_pInt mesh_element(10,e) = mesh_element(9,e) + 1_pInt
mesh_element(11,e) = mesh_element(9,e) + res(1) + 1_pInt mesh_element(11,e) = mesh_element(9,e) + res(1) + 2_pInt
mesh_element(12,e) = mesh_element(9,e) + res(1) mesh_element(12,e) = mesh_element(9,e) + res(1) + 1_pInt
mesh_maxValStateVar(1) = max(mesh_maxValStateVar(1),mesh_element(3,e)) !needed for statistics mesh_maxValStateVar(1) = max(mesh_maxValStateVar(1),mesh_element(3,e)) !needed for statistics
mesh_maxValStateVar(2) = max(mesh_maxValStateVar(2),mesh_element(4,e)) mesh_maxValStateVar(2) = max(mesh_maxValStateVar(2),mesh_element(4,e))
enddo enddo
@ -3850,7 +3845,7 @@ function mesh_spectral_getDimension(fileUnit)
if (keyword(1:4) == 'head') then if (keyword(1:4) == 'head') then
headerLength = IO_intValue(line,positions,1_pInt) + 1_pInt headerLength = IO_intValue(line,positions,1_pInt) + 1_pInt
else else
call IO_error(error_ID=842_pInt) call IO_error(error_ID=841_pInt, ext_msg='mesh_spectral_getDimension')
endif endif
rewind(myUnit) rewind(myUnit)
do i = 1_pInt, headerLength do i = 1_pInt, headerLength
@ -3877,7 +3872,7 @@ function mesh_spectral_getDimension(fileUnit)
if (.not. gotDimension) & if (.not. gotDimension) &
call IO_error(error_ID = 845_pInt, ext_msg='dimension') call IO_error(error_ID = 845_pInt, ext_msg='dimension')
if (any(mesh_spectral_getDimension<=0.0_pReal)) & if (any(mesh_spectral_getDimension<=0.0_pReal)) &
call IO_error(error_ID = 802_pInt, ext_msg='dimension') call IO_error(error_ID = 844_pInt, ext_msg='mesh_spectral_getDimension')
end function mesh_spectral_getDimension end function mesh_spectral_getDimension
@ -3922,7 +3917,7 @@ function mesh_spectral_getResolution(fileUnit)
if (keyword(1:4) == 'head') then if (keyword(1:4) == 'head') then
headerLength = IO_intValue(line,positions,1_pInt) + 1_pInt headerLength = IO_intValue(line,positions,1_pInt) + 1_pInt
else else
call IO_error(error_ID=842_pInt) call IO_error(error_ID=841_pInt, ext_msg='mesh_spectral_getResolution')
endif endif
rewind(myUnit) rewind(myUnit)
do i = 1_pInt, headerLength do i = 1_pInt, headerLength
@ -3948,11 +3943,14 @@ function mesh_spectral_getResolution(fileUnit)
if (.not. gotResolution) & if (.not. gotResolution) &
call IO_error(error_ID = 845_pInt, ext_msg='resolution') call IO_error(error_ID = 845_pInt, ext_msg='resolution')
if(mod(mesh_spectral_getResolution(1),2_pInt)/=0_pInt .or.& if((mod(mesh_spectral_getResolution(1),2_pInt)/=0_pInt .or. & ! must be a even number
mod(mesh_spectral_getResolution(2),2_pInt)/=0_pInt .or.& mesh_spectral_getResolution(1) < 2_pInt .or. & ! and larger than 1
mod(mesh_spectral_getResolution(2),2_pInt)/=0_pInt .or. & ! -"-
mesh_spectral_getResolution(2) < 2_pInt .or. & ! -"-
(mod(mesh_spectral_getResolution(3),2_pInt)/=0_pInt .and. & (mod(mesh_spectral_getResolution(3),2_pInt)/=0_pInt .and. &
mesh_spectral_getResolution(3)/= 1_pInt))& mesh_spectral_getResolution(3)/= 1_pInt)) .or. & ! third res might be 1
call IO_error(error_ID = 802_pInt, ext_msg='resolution') mesh_spectral_getResolution(3) < 1_pInt) &
call IO_error(error_ID = 843_pInt, ext_msg='mesh_spectral_getResolution')
end function mesh_spectral_getResolution end function mesh_spectral_getResolution
@ -3995,7 +3993,7 @@ function mesh_spectral_getHomogenization(fileUnit)
if (keyword(1:4) == 'head') then if (keyword(1:4) == 'head') then
headerLength = IO_intValue(line,positions,1_pInt) + 1_pInt headerLength = IO_intValue(line,positions,1_pInt) + 1_pInt
else else
call IO_error(error_ID=842_pInt) call IO_error(error_ID=841_pInt, ext_msg='mesh_spectral_getHomogenization')
endif endif
rewind(myUnit) rewind(myUnit)
do i = 1_pInt, headerLength do i = 1_pInt, headerLength
@ -4013,7 +4011,7 @@ function mesh_spectral_getHomogenization(fileUnit)
if (.not. gotHomogenization ) & if (.not. gotHomogenization ) &
call IO_error(error_ID = 845_pInt, ext_msg='homogenization') call IO_error(error_ID = 845_pInt, ext_msg='homogenization')
if (mesh_spectral_getHomogenization<1_pInt) & if (mesh_spectral_getHomogenization<1_pInt) &
call IO_error(error_ID = 802_pInt, ext_msg='homogenization') call IO_error(error_ID = 842_pInt, ext_msg='mesh_spectral_getHomogenization')
end function mesh_spectral_getHomogenization end function mesh_spectral_getHomogenization