From c752dd5474f3f34f264768c2808b10510d3426e2 Mon Sep 17 00:00:00 2001 From: Krishna Komerla Date: Tue, 8 May 2012 14:57:06 +0000 Subject: [PATCH] regridding is now working, changed the subroutine into a function changed order of arrays in nearest neighbor search to make it fortran fast constitutive.f90 and homogenization.f90 write state size out during initialization setup/setup_processing.py is using byterecl to be compatible with binary files written out by solver --- code/DAMASK_spectral_interface.f90 | 12 +- code/constitutive.f90 | 8 + code/damask.core.pyf | 11 +- code/homogenization.f90 | 11 +- code/math.f90 | 16 +- code/mesh.f90 | 258 ++++++++++++++++----------- processing/setup/setup_processing.py | 6 +- 7 files changed, 194 insertions(+), 128 deletions(-) diff --git a/code/DAMASK_spectral_interface.f90 b/code/DAMASK_spectral_interface.f90 index 0582e14bb..f0337447c 100644 --- a/code/DAMASK_spectral_interface.f90 +++ b/code/DAMASK_spectral_interface.f90 @@ -183,12 +183,12 @@ subroutine DAMASK_interface_init(loadcaseParameterIn,geometryParameterIn) call GET_ENVIRONMENT_VARIABLE('HOST',hostName) call GET_ENVIRONMENT_VARIABLE('USER',userName) - write(6,'(a,2(i2.2,a),i4.4)') ' Date: ',dateAndTime(3),'/',& - dateAndTime(2),'/',& - dateAndTime(1) - write(6,'(a,2(i2.2,a),i2.2)') ' Time: ',dateAndTime(5),':',& - dateAndTime(6),':',& - dateAndTime(7) + write(6,'(a,2(i2.2,a),i4.4)') 'Date: ',dateAndTime(3),'/',& + dateAndTime(2),'/',& + dateAndTime(1) + write(6,'(a,2(i2.2,a),i2.2)') 'Time: ',dateAndTime(5),':',& + dateAndTime(6),':',& + dateAndTime(7) write(6,'(a,a)') 'Host Name: ', trim(hostName) write(6,'(a,a)') 'User Name: ', trim(userName) write(6,'(a,a)') 'Path Separator: ', getPathSep() diff --git a/code/constitutive.f90 b/code/constitutive.f90 index aa87f9ae4..49aef7958 100644 --- a/code/constitutive.f90 +++ b/code/constitutive.f90 @@ -29,6 +29,7 @@ MODULE constitutive use prec, only: pInt, p_vec +use IO, only: IO_write_jobBinaryFile implicit none type(p_vec), dimension(:,:,:), allocatable :: & @@ -356,6 +357,13 @@ endif enddo !$OMP END PARALLEL DO +!----- write out state size file---------------- +open(777) +call IO_write_jobBinaryFile(777,'sizeStateConst', size(constitutive_sizeState)) +write (777,rec=1) constitutive_sizeState +close(777) +!----------------------------------------------- + constitutive_maxSizeState = maxval(constitutive_sizeState) constitutive_maxSizeDotState = maxval(constitutive_sizeDotState) constitutive_maxSizePostResults = maxval(constitutive_sizePostResults) diff --git a/code/damask.core.pyf b/code/damask.core.pyf index a439df716..7f6156e3c 100644 --- a/code/damask.core.pyf +++ b/code/damask.core.pyf @@ -167,8 +167,8 @@ python module core ! in real*8, dimension(3), intent(in) :: geomdim integer, intent(in) :: queryPoints integer, intent(in) :: domainPoints - real*8, dimension(queryPoints,spatialDim), intent(in), depend(queryPoints,spatialDim) :: querySet - real*8, dimension(domainPoints,spatialDim), intent(in), depend(domainPoints,spatialDim) :: domainSet + real*8, dimension(spatialDim,queryPoints), intent(in), depend(queryPoints,spatialDim) :: querySet + real*8, dimension(spatialDim,domainPoints), intent(in), depend(domainPoints,spatialDim) :: domainSet ! output variable integer, dimension(queryPoints), intent(out), depend(queryPoints) :: indices ! depending on input @@ -178,9 +178,10 @@ python module core ! in module mesh ! in :mesh:mesh.f90 - subroutine mesh_regrid(resNew) ! in :mesh:mesh.f90 - integer, dimension(3), intent(in,out), optional :: resNew - end subroutine mesh_regrid + function mesh_regrid(resNewInput) ! in :mesh:mesh.f90 + integer, dimension(3) :: mesh_regrid + integer, dimension(3), intent(in), optional :: resNewInput + end function mesh_regrid end module mesh diff --git a/code/homogenization.f90 b/code/homogenization.f90 index 778cee9fd..ed327669c 100644 --- a/code/homogenization.f90 +++ b/code/homogenization.f90 @@ -34,6 +34,7 @@ MODULE homogenization !*** Include other modules *** use prec, only: pInt,pReal,p_vec + use IO, only: IO_write_jobBinaryFile implicit none ! **************************************************************** @@ -161,7 +162,7 @@ forall (i = 1:mesh_maxNips,e = 1:mesh_NcpElems) end forall -! --- ALLOCATE AND INITIALIZE GLOBAL STATE AND POSTRESULTS VARIABLES --- +! --- ALLOCATE AND INITIALIZE GLOBAL STATE AND POSTRESULTS VARIABLES ---------- !$OMP PARALLEL DO PRIVATE(myInstance) do e = 1,mesh_NcpElems ! loop over elements @@ -192,6 +193,14 @@ end forall enddo enddo !$OMP END PARALLEL DO + +!---write state size file out--------------------------------------- +open(777) +call IO_write_jobBinaryFile(777,'sizeStateHomog',size(homogenization_sizeState)) +write (777,rec=1) homogenization_sizeState +close(777) +!-------------------------------------------------------------- + homogenization_maxSizeState = maxval(homogenization_sizeState) homogenization_maxSizePostResults = maxval(homogenization_sizePostResults) materialpoint_sizeResults = 1 & ! grain count diff --git a/code/math.f90 b/code/math.f90 index ab7ff6f7b..11a5b45d0 100644 --- a/code/math.f90 +++ b/code/math.f90 @@ -3832,19 +3832,19 @@ subroutine math_nearestNeighborSearch(spatialDim, Favg, geomdim, queryPoints, do real(pReal), dimension(3), intent(in) :: geomdim integer(pInt), intent(in) :: domainPoints integer(pInt), intent(in) :: queryPoints - real(pReal), dimension(queryPoints,spatialDim), intent(in) :: querySet - real(pReal), dimension(domainPoints,spatialDim), intent(in) :: domainSet + real(pReal), dimension(spatialDim,queryPoints), intent(in) :: querySet + real(pReal), dimension(spatialDim,domainPoints), intent(in) :: domainSet ! output variable integer(pInt), dimension(queryPoints), intent(out) :: indices ! other variables depending on input - real(pReal), dimension(3,(3_pInt**spatialDim)*domainPoints) :: domainSetLarge + real(pReal), dimension(spatialDim,(3_pInt**spatialDim)*domainPoints) :: domainSetLarge ! other variables integer(pInt) :: i,j, l,m,n type(kdtree2), pointer :: tree type(kdtree2_result), dimension(1) :: Results - if (size(querySet(1,:)) /= spatialDim) call IO_error(407_pInt,ext_msg='query set') - if (size(domainSet(1,:)) /= spatialDim) call IO_error(407_pInt,ext_msg='domain set') + if (size(querySet(:,1)) /= spatialDim) call IO_error(407_pInt,ext_msg='query set') + if (size(domainSet(:,1)) /= spatialDim) call IO_error(407_pInt,ext_msg='domain set') i = 0_pInt @@ -3852,14 +3852,14 @@ subroutine math_nearestNeighborSearch(spatialDim, Favg, geomdim, queryPoints, do do j = 1_pInt, domainPoints do l = -1_pInt, 1_pInt; do m = -1_pInt, 1_pInt i = i + 1_pInt - domainSetLarge(1:3,i) = domainSet(j,1:3) + math_mul33x3(Favg,real([l,m,0_pInt],pReal)*geomdim) + domainSetLarge(1:2,i) = domainSet(1:2,j) +matmul(Favg(1:2,1:2),real([l,m],pReal)*geomdim(1:2)) enddo; enddo enddo else do j = 1_pInt, domainPoints do l = -1_pInt, 1_pInt; do m = -1_pInt, 1_pInt; do n = -1_pInt, 1_pInt i = i + 1_pInt - domainSetLarge(1:3,i) = domainSet(j,1:3) + math_mul33x3(Favg,real([l,m,n],pReal)*geomdim) + domainSetLarge(1:3,i) = domainSet(1:3,j) + math_mul33x3(Favg,real([l,m,n],pReal)*geomdim) enddo; enddo; enddo enddo endif @@ -3867,7 +3867,7 @@ subroutine math_nearestNeighborSearch(spatialDim, Favg, geomdim, queryPoints, do tree => kdtree2_create(domainSetLarge,sort=.true.,rearrange=.true.) do j = 1_pInt, queryPoints - call kdtree2_n_nearest(tp=tree, qv=querySet(j,1:spatialDim),nn=1_pInt, results = Results) + call kdtree2_n_nearest(tp=tree, qv=querySet(1:spatialDim,j),nn=1_pInt, results = Results) indices(j) = Results(1)%idx enddo diff --git a/code/mesh.f90 b/code/mesh.f90 index b1709dd40..e279d6918 100644 --- a/code/mesh.f90 +++ b/code/mesh.f90 @@ -3528,7 +3528,7 @@ end subroutine mesh_tell_statistics -subroutine mesh_regrid(resNew) !use new_res=[0.0,0.0,0.0] for automatic determination of new grid +function mesh_regrid(resNewInput) !use new_res=[0.0,0.0,0.0] for automatic determination of new grid use prec, only: & pInt, & pReal @@ -3542,86 +3542,122 @@ subroutine mesh_regrid(resNew) !use new_res=[0.0,0.0,0.0] for automati deformed_FFT, & math_mul33x3 - integer(pInt), dimension(3), optional, intent(inout) :: resNew - integer(pInt):: maxsize, i, j, k, ielem, Npoints, NpointsNew - integer(pInt), dimension(3) :: res - integer(pInt), dimension(:), allocatable :: indices, outputSize - real(pReal), dimension(3) :: geomdim = 0.0_pReal + integer(pInt), dimension(3) :: mesh_regrid + integer(pInt), dimension(3), optional, intent(in) :: resNewInput + integer(pInt):: maxsize, i, j, k, ielem, Npoints, NpointsNew, spatialDim + integer(pInt), dimension(3) :: res, resNew + integer(pInt), dimension(:), allocatable :: indices + real(pReal), dimension(3) :: geomdim real(pReal), dimension(3,3) :: Favg real(pReal), dimension(:,:), allocatable :: & - coordinatesNew + coordinatesNew, & + coordinatesLinear + + real(pReal), dimension(:,:,:), allocatable :: & + material_phase, material_phaseNew + real(pReal), dimension(:,:,:,:,:), allocatable :: & F, FNew, & Fp, FpNew, & Lp, LpNew, & dcsdE, dcsdENew + real(pReal), dimension (:,:,:,:,:,:,:), allocatable :: & dPdF, dPdFNew real(pReal), dimension (:,:,:,:), allocatable :: & - coordinates, & - Tstar, TstarNew, & - stateConst, stateConstNew, & - stateHomog, stateHomogNew - integer(pInt), dimension (:,:), allocatable :: SizeConst + coordinates, & + Tstar, TstarNew, & + stateHomog, & + stateConst + integer(pInt), dimension(:,:), allocatable :: & + sizeStateConst, sizeStateHomog + res = mesh_spectral_getResolution() geomdim = mesh_spectral_getDimension() - - - if (.not. present(resNew)) then - resNew=res - endif - Npoints = res(1)*res(2)*res(3) - NpointsNew = resNew(1)*resNew(2)*resNew(3) - print*, 'resolution ', res - print*, 'new resolution', resNew + - allocate(coordinates(res(1),res(2),res(3),3)) - allocate(coordinatesNew(resNew(1)*resNew(2)*resNew(3),3)) - allocate(indices(Npoints)) allocate(F(res(1),res(2),res(3),3,3)) - - call IO_read_jobBinaryFile(777,'convergedSpectralDefgrad',trim(getSolverJobName()),size(F)) read (777,rec=1) F close (777) -! ----Calculate average deformation for remesh-------- +! ----Calculate deformed configuration and average-------- do i= 1_pInt,3_pInt; do j = 1_pInt,3_pInt Favg(i,j) = sum(F(1:res(1),1:res(2),1:res(3),i,j)) / Npoints enddo; enddo - + allocate(coordinates(res(1),res(2),res(3),3)) call deformed_fft(res,geomdim,Favg,1.0_pReal,F,coordinates) + deallocate(F) +! ----Store coordinates into a linear list-------- + allocate(coordinatesLinear(3,Npoints)) + ielem = 0_pInt + do k=1_pInt,res(3); do j=1_pInt, res(2); do i=1_pInt, res(1) + ielem = ielem + 1_pInt + coordinatesLinear(1:3,ielem) = coordinates(i,j,k,1:3) + enddo; enddo; enddo + deallocate(coordinates) + + ! ----For 2D /3D case---------------------------------- + if (res(3)== 1_pInt) then + spatialDim = 2_pInt + else + spatialDim = 3_pInt + endif + +! ----determine/calculate new resolution-------------- + if (.not. present(resNewInput)) then + resNew=res + elseif(all(resNewInput==0.0_pReal)) then + !do automatic determination + else + resNew=resNewInput + endif + NpointsNew = resNew(1)*resNew(2)*resNew(3) + +! ----Calculate regular new coordinates------------ + allocate(coordinatesNew(3,NpointsNew)) ielem = 0_pInt do k=1_pInt,resNew(3); do j=1_pInt, resNew(2); do i=1_pInt, resNew(1) ielem = ielem + 1_pInt - coordinatesNew(ielem,1:3) = math_mul33x3(Favg, geomdim/real(resNew,pReal)*real([i,j,k],pReal) & + coordinatesNew(1:3,ielem) = math_mul33x3(Favg, geomdim/real(resNew,pReal)*real([i,j,k],pReal) & - geomdim/real(2_pInt*resNew,pReal)) enddo; enddo; enddo !----- Nearest neighbour search ----------------------- - call math_nearestNeighborSearch(res, Favg, geomdim, NpointsNew, & - coordinates, coordinatesNew, indices) - indices = indices /3_pInt**3_pInt +1_pInt ! adjust it for 2D case + allocate(indices(Npoints)) + call math_nearestNeighborSearch(spatialDim, Favg, geomdim, NpointsNew, Npoints, & + coordinatesNew, coordinatesLinear, indices) + deallocate(coordinatesNew) + indices = indices / 3_pInt**spatialDim +1_pInt - deallocate(F) - deallocate(coordinates) - deallocate(coordinatesNew) +!--------------------------------------------------------------------------- + allocate(material_phase (1,1, Npoints)) + allocate(material_phaseNew (1,1, NpointsNew)) + call IO_read_jobBinaryFile(777,'recordedPhase',trim(getSolverJobName()),size(material_phase)) + read (777,rec=1) material_phase + close (777) + do i = 1, NpointsNew + material_phaseNew(1,1,i) = material_phase(1,1,indices(i)) + enddo -! --------------------------------------------------------------------------- - allocate(F (3,3,1,1, Npoints)) - allocate(FNew(3,3,1,1, NpointsNew)) + call IO_write_jobBinaryFile(777,'recordedPhase',size(material_phaseNew)) + write (777,rec=1) material_phaseNew + close (777) + deallocate(material_phase) + deallocate(material_phaseNew) + +!--------------------------------------------------------------------------- + allocate(F (3,3,1,1, Npoints)) + allocate(FNew (3,3,1,1, NpointsNew)) call IO_read_jobBinaryFile(777,'convergedF',trim(getSolverJobName()),size(F)) read (777,rec=1) F close (777) - call flush(6) - do i = 1, NpointsNew FNew(1:3,1:3,1,1,i) = F(1:3,1:3,1,1,indices(i)) enddo - call flush(6) call IO_write_jobBinaryFile(777,'convergedF',size(FNew)) write (777,rec=1) FNew @@ -3630,15 +3666,14 @@ subroutine mesh_regrid(resNew) !use new_res=[0.0,0.0,0.0] for automati deallocate(FNew) !--------------------------------------------------------------------- - allocate(Fp (3,3,1,1,Npoints)) - allocate(FpNew(3,3,1,1,NpointsNew)) + allocate(Fp (3,3,1,1,Npoints)) + allocate(FpNew (3,3,1,1,NpointsNew)) call IO_read_jobBinaryFile(777,'convergedFp',trim(getSolverJobName()),size(Fp)) read (777,rec=1) Fp close (777) - do i = 1, NpointsNew FpNew(1:3,1:3,1,1,i) = Fp(1:3,1:3,1,1,indices(i)) - enddo + enddo call IO_write_jobBinaryFile(777,'convergedFp',size(FpNew)) write (777,rec=1) FpNew @@ -3647,50 +3682,44 @@ subroutine mesh_regrid(resNew) !use new_res=[0.0,0.0,0.0] for automati deallocate(FpNew) !------------------------------------------------------------------------ - allocate(Lp (3,3,1,1,Npoints)) - allocate(LpNew(3,3,1,1,NpointsNew)) + allocate(Lp (3,3,1,1,Npoints)) + allocate(LpNew (3,3,1,1,NpointsNew)) call IO_read_jobBinaryFile(777,'convergedLp',trim(getSolverJobName()),size(Lp)) read (777,rec=1) Lp close (777) - do i = 1, NpointsNew LpNew(1:3,1:3,1,1,i) = Lp(1:3,1:3,1,1,indices(i)) enddo - call IO_write_jobBinaryFile(777,'convergedLp',size(LpNew)) write (777,rec=1) LpNew close (777) deallocate(Lp) deallocate(LpNew) - !---------------------------------------------------------------------------- - allocate(dcsdE (6,6,1,1,Npoints)) - allocate(dcsdENew(6,6,1,1,NpointsNew)) +!---------------------------------------------------------------------------- + allocate(dcsdE (6,6,1,1,Npoints)) + allocate(dcsdENew (6,6,1,1,NpointsNew)) call IO_read_jobBinaryFile(777,'convergeddcsdE',trim(getSolverJobName()),size(dcsdE)) read (777,rec=1) dcsdE close (777) - do i = 1, NpointsNew dcsdENew(1:6,1:6,1,1,i) = dcsdE(1:6,1:6,1,1,indices(i)) enddo - call IO_write_jobBinaryFile(777,'convergeddcsdE',size(dcsdENew)) write (777,rec=1) dcsdENew close (777) deallocate(dcsdE) deallocate(dcsdENew) -! --------------------------------------------------------------------------- - allocate(dPdF (3,3,3,3,1,1,Npoints)) - allocate(dPdFNew(3,3,3,3,1,1,NpointsNew)) +!--------------------------------------------------------------------------- + allocate(dPdF (3,3,3,3,1,1,Npoints)) + allocate(dPdFNew (3,3,3,3,1,1,NpointsNew)) call IO_read_jobBinaryFile(777,'convergeddPdF',trim(getSolverJobName()),size(dPdF)) read (777,rec=1) dPdF close (777) - do i = 1, NpointsNew dPdFNew(1:3,1:3,1:3,1:3,1,1,i) = dPdF(1:3,1:3,1:3,1:3,1,1,indices(i)) enddo - call IO_write_jobBinaryFile(777,'convergeddPdF',size(dPdFNew)) write (777,rec=1) dPdFNew close (777) @@ -3698,66 +3727,85 @@ subroutine mesh_regrid(resNew) !use new_res=[0.0,0.0,0.0] for automati deallocate(dPdFNew) !--------------------------------------------------------------------------- - allocate(Tstar (6,1,1,Npoints)) - allocate(TstarNew(6,1,1,NpointsNew)) + allocate(Tstar (6,1,1,Npoints)) + allocate(TstarNew (6,1,1,NpointsNew)) call IO_read_jobBinaryFile(777,'convergedTstar',trim(getSolverJobName()),size(Tstar)) read (777,rec=1) Tstar close (777) - do i = 1, NpointsNew TstarNew(1:6,1,1,i) = Tstar(1:6,1,1,indices(i)) enddo - call IO_write_jobBinaryFile(777,'convergedTstar',size(TstarNew)) write (777,rec=1) TstarNew close (777) deallocate(Tstar) deallocate(TstarNew) + !---------------------------------------------------------------------------- -!for homog and state const - ! allocate(SizeConst (1,Npoints)) - ! allocate(StateConst (1,Npoints)) - ! allocate(StateConstNew (1,Npoints)) - ! call IO_read_jobBinaryFile(777,'SizeConst',trim(getSolverJobName())) - ! read (777,rec=1) SizeConst - ! close(777) - ! maxsize = maxval(SizeConst) - ! call IO_read_jobBinaryFile(777,'convergedStateConst',trim(getSolverJobName())) - ! do i =1, Npoints - ! do l = 1,SizeConst(1,i) - ! read(777,rec=i) StateConst(1,1,i,l) - ! enddo - ! enddo - ! close(777) + allocate(sizeStateConst(1,Npoints)) + call IO_read_jobBinaryFile(777,'sizeStateConst',trim(getSolverJobName()),size(sizeStateConst)) + read (777,rec=1) sizeStateConst + close (777) + maxsize = maxval(sizeStateConst) + allocate(StateConst (1,1,Npoints,maxsize)) + + call IO_read_jobBinaryFile(777,'convergedStateConst',trim(getSolverJobName())) + k = 0_pInt + do i =1, Npoints + do j = 1,sizeStateConst(1,i) + k = k+1_pInt + read(777,rec=k) StateConst(1,1,i,j) + enddo + enddo + close(777) + + call IO_write_jobBinaryFile(777,'convergedStateConst') + k = 0_pInt + do i = 1,NpointsNew + do j = 1,sizeStateConst(1,i) + k=k+1_pInt + write(777,rec=k) StateConst(1,1,indices(i),j) + enddo + enddo + close (777) + deallocate(sizeStateConst) + deallocate(StateConst) + +!---------------------------------------------------------------------------- + allocate(sizeStateHomog(1,Npoints)) + call IO_read_jobBinaryFile(777,'sizeStateHomog',trim(getSolverJobName()),size(sizeStateHomog)) + read (777,rec=1) sizeStateHomog + close (777) + maxsize = maxval(sizeStateHomog) + allocate(stateHomog (1,1,Npoints,maxsize)) + + call IO_read_jobBinaryFile(777,'convergedStateHomog',trim(getSolverJobName())) + k = 0_pInt + do i =1, Npoints + do j = 1,sizeStateHomog(1,i) + k = k+1_pInt + read(777,rec=k) stateHomog(1,1,i,j) + enddo + enddo + close(777) + + call IO_write_jobBinaryFile(777,'convergedStateHomog') + k = 0_pInt + do i = 1,NpointsNew + do j = 1,sizeStateHomog(1,i) + k=k+1_pInt + write(777,rec=k) stateHomog(1,1,indices(i),j) + enddo + enddo + close (777) + deallocate(sizeStateHomog) + deallocate(stateHomog) - ! do i = 1, NpointsNew - ! do l = 1,SizeConst(1,i) - ! StateConstNew(1,1,i,l) = StateConst(1,1,indices(i),l) - ! enddo - ! enddo + deallocate(indices) - ! call IO_write_jobBinaryFile(777,'convergedStateConst',trim(getSolverJobName())) - ! m = 0_pInt - ! do i = 1,NpointsNew - ! do l = 1,SizeConst(1,i) - ! write(777,rec=m) StateConstNew(1,1,i,l) - ! enddo - ! enddo - ! close (777) - !---------------------------------------------------------------------------- - ! State homog - ! call IO_read_jobBinaryFile(777,'convergedStateHomog',FEmodelGeometry) - ! m = 0_pInt - ! do k = 1,Npoints; - ! do l = 1,homogenization_sizeState(j,k) - ! m = m+1_pInt - ! read(777,rec=m) stateHomog(1,k)%p(l) - ! enddo - ! enddo; enddo - ! close (777) - - -end subroutine mesh_regrid + mesh_regrid = resNew +end function mesh_regrid + function mesh_spectral_getDimension(fileUnit) use IO, only: & diff --git a/processing/setup/setup_processing.py b/processing/setup/setup_processing.py index 6f89afc1a..6ae6f4f8d 100755 --- a/processing/setup/setup_processing.py +++ b/processing/setup/setup_processing.py @@ -50,9 +50,9 @@ if options.compiler not in compilers: f2py_compiler = { 'gfortran': 'gnu95 --f90flags="-fno-range-check -xf95-cpp-input -std=f2008 -fall-intrinsics"', 'gnu95': 'gnu95 --f90flags="-fno-range-check -xf95-cpp-input -std=f2008 -fall-intrinsics"', - 'intel32': 'intel --f90flags="-fpp -stand f03 -diag-disable 5268"', - 'intel': 'intelem --f90flags="-fpp -stand f03 -diag-disable 5268"', - 'ifort': 'intelem --f90flags="-fpp -stand f03 -diag-disable 5268"', + 'intel32': 'intel --f90flags="-fpp -stand f03 -diag-disable 5268 -assume byterecl"', + 'intel': 'intelem --f90flags="-fpp -stand f03 -diag-disable 5268 -assume byterecl"', + 'ifort': 'intelem --f90flags="-fpp -stand f03 -diag-disable 5268 -assume byterecl"', }[options.compiler] compiler = { 'gfortran': 'gfortran',