moved calculation and output of geometry related data to mesh.f90, moved output of command line related information to DAMASK_spectral_interface.f90.
improved new solver structure (use with make SOLVER=NEW)
This commit is contained in:
Martin Diehl
parent
a48d180f8c
commit
7e3a837640
|
|
@ -77,10 +77,14 @@ program DAMASK_spectral
|
|||
use math
|
||||
|
||||
use mesh, only : &
|
||||
mesh_spectral_getResolution, &
|
||||
mesh_spectral_getDimension, &
|
||||
mesh_spectral_getHomogenization
|
||||
|
||||
homog, &
|
||||
res, &
|
||||
res1_red, &
|
||||
mesh_NcpElems, &
|
||||
wgt, &
|
||||
geomdim, &
|
||||
virt_dim
|
||||
|
||||
use CPFEM, only: &
|
||||
CPFEM_general, &
|
||||
CPFEM_initAll
|
||||
|
|
@ -96,8 +100,7 @@ program DAMASK_spectral
|
|||
rotation_tol, &
|
||||
itmax,&
|
||||
itmin, &
|
||||
memory_efficient, &
|
||||
divergence_correction, &
|
||||
memory_efficient, &
|
||||
DAMASK_NumThreadsInt, &
|
||||
fftw_planner_flag, &
|
||||
fftw_timelimit
|
||||
|
|
@ -132,10 +135,7 @@ program DAMASK_spectral
|
|||
N_l = 0_pInt, &
|
||||
N_t = 0_pInt, &
|
||||
N_n = 0_pInt, &
|
||||
N_Fdot = 0_pInt, &
|
||||
Npoints,& ! number of Fourier points
|
||||
homog, & ! homogenization scheme used
|
||||
res1_red ! to store res(1)/2 +1
|
||||
N_Fdot = 0_pInt
|
||||
|
||||
character(len=1024) :: &
|
||||
line
|
||||
|
|
@ -159,11 +159,6 @@ program DAMASK_spectral
|
|||
|
||||
type(bc_type), allocatable, dimension(:) :: bc
|
||||
|
||||
|
||||
real(pReal) :: wgt
|
||||
real(pReal), dimension(3) :: geomdim = 0.0_pReal, virt_dim = 0.0_pReal ! physical dimension of volume element per direction
|
||||
integer(pInt), dimension(3) :: res = 1_pInt ! resolution (number of Fourier points) in each direction
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! stress, stiffness and compliance average etc.
|
||||
real(pReal), dimension(3,3) :: &
|
||||
|
|
@ -368,15 +363,6 @@ program DAMASK_spectral
|
|||
end select
|
||||
enddo; enddo
|
||||
101 close(myUnit)
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! get resolution, dimension, homogenization and variables derived from resolution
|
||||
res = mesh_spectral_getResolution()
|
||||
geomdim = mesh_spectral_getDimension()
|
||||
homog = mesh_spectral_getHomogenization()
|
||||
res1_red = res(1)/2_pInt + 1_pInt ! size of complex array in first dimension (c2r, r2c)
|
||||
Npoints = res(1)*res(2)*res(3)
|
||||
wgt = 1.0_pReal/real(Npoints, pReal)
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! output of geometry
|
||||
|
|
@ -541,8 +527,8 @@ program DAMASK_spectral
|
|||
call IO_read_jobBinaryFile(777,'F_aim_lastInc',trim(getSolverJobName()),size(F_aim_lastInc))
|
||||
read (777,rec=1) F_aim_lastInc
|
||||
close (777)
|
||||
|
||||
coordinates = 0.0 ! change it later!!!
|
||||
CPFEM_mode = 2_pInt
|
||||
endif
|
||||
ielem = 0_pInt
|
||||
do k = 1_pInt, res(3); do j = 1_pInt, res(2); do i = 1_pInt, res(1)
|
||||
|
|
@ -561,13 +547,6 @@ program DAMASK_spectral
|
|||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! calculation of discrete angular frequencies, ordered as in FFTW (wrap around)
|
||||
if (divergence_correction) then
|
||||
do i = 1_pInt, 3_pInt
|
||||
if (i /= minloc(geomdim,1) .and. i /= maxloc(geomdim,1)) virt_dim = geomdim/geomdim(i)
|
||||
enddo
|
||||
else
|
||||
virt_dim = geomdim
|
||||
endif
|
||||
|
||||
do k = 1_pInt, res(3)
|
||||
k_s(3) = k - 1_pInt
|
||||
|
|
@ -632,7 +611,7 @@ program DAMASK_spectral
|
|||
write(538) 'increments', bc(1:N_Loadcases)%incs ! one entry per loadcase
|
||||
write(538) 'startingIncrement', restartInc - 1_pInt ! start with writing out the previous inc
|
||||
write(538) 'eoh' ! end of header
|
||||
write(538) materialpoint_results(1_pInt:materialpoint_sizeResults,1,1_pInt:Npoints) ! initial (non-deformed or read-in) results
|
||||
write(538) materialpoint_results(1_pInt:materialpoint_sizeResults,1,1_pInt:mesh_NcpElems) ! initial (non-deformed or read-in) results
|
||||
if (debugGeneral) write(6,'(a)') 'Header of result file written out'
|
||||
endif
|
||||
|
||||
|
|
@ -661,7 +640,7 @@ program DAMASK_spectral
|
|||
!##################################################################################################
|
||||
do inc = 1_pInt, bc(loadcase)%incs
|
||||
totalIncsCounter = totalIncsCounter + 1_pInt
|
||||
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! forwarding time
|
||||
timeinc_old = timeinc
|
||||
|
|
@ -711,7 +690,7 @@ program DAMASK_spectral
|
|||
F_lastInc(i,j,k,1:3,1:3) = temp33_Real
|
||||
Favg = Favg + F(i,j,k,1:3,1:3)
|
||||
enddo; enddo; enddo
|
||||
deltaF_aim = guessmode *(Favg*wgt -F_aim) ! average correction in case of guessing to
|
||||
deltaF_aim = guessmode *(Favg*wgt -math_rotate_backward33(F_aim,bc(loadcase)%rotation)) ! average correction in case of guessing to
|
||||
do k = 1_pInt, res(3); do j = 1_pInt, res(2); do i = 1_pInt, res(1)
|
||||
F(i,j,k,1:3,1:3) = F(i,j,k,1:3,1:3) - deltaF_aim ! correct in case avg of F is not F_aim
|
||||
enddo; enddo; enddo
|
||||
|
|
@ -806,9 +785,9 @@ program DAMASK_spectral
|
|||
restartWrite = .false.
|
||||
|
||||
! for test of regridding
|
||||
! if( bc(loadcase)%restartFrequency > 0_pInt .and. &
|
||||
! mod(inc-1,bc(loadcase)%restartFrequency) == 0_pInt .and. &
|
||||
! restartInc/=inc) call quit(-1*(restartInc+1)) ! trigger exit to regrid
|
||||
! if( bc(loadcase)%restartFrequency > 0_pInt .and. &
|
||||
! mod(inc-1,bc(loadcase)%restartFrequency) == 0_pInt .and. &
|
||||
! restartInc/=inc) call quit(-1*(restartInc+1)) ! trigger exit to regrid
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! copy one component of the stress field to to a single FT and check for mismatch
|
||||
|
|
@ -981,7 +960,7 @@ program DAMASK_spectral
|
|||
|
||||
endif
|
||||
deltaF_fourier(1,1,1,1:3,1:3) = cmplx((F_aim_lab_lastIter - F_aim_lab) & ! assign (negative) average deformation gradient change to zero frequency (real part)
|
||||
* real(Npoints,pReal),0.0_pReal,pReal) ! singular point at xi=(0.0,0.0,0.0) i.e. i=j=k=1
|
||||
* real(mesh_NcpElems,pReal),0.0_pReal,pReal) ! singular point at xi=(0.0,0.0,0.0) i.e. i=j=k=1
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! comparing 1 and 3x3 inverse FT results
|
||||
|
|
@ -1077,7 +1056,7 @@ program DAMASK_spectral
|
|||
if (mod(inc,bc(loadcase)%outputFrequency) == 0_pInt) then ! at output frequency
|
||||
write(6,'(a)') ''
|
||||
write(6,'(a)') '... writing results to file ......................................'
|
||||
write(538) materialpoint_results(1_pInt:materialpoint_sizeResults,1,1_pInt:Npoints) ! write result to file
|
||||
write(538) materialpoint_results(1_pInt:materialpoint_sizeResults,1,1_pInt:mesh_NcpElems) ! write result to file
|
||||
flush(538)
|
||||
endif
|
||||
|
||||
|
|
@ -1092,9 +1071,6 @@ program DAMASK_spectral
|
|||
call IO_write_jobBinaryFile(777,'convergedSpectralDefgrad_lastInc',size(F_lastInc)) ! writing F_lastInc field to file
|
||||
write (777,rec=1) F_lastInc
|
||||
close (777)
|
||||
call IO_write_jobBinaryFile(777,'C',size(C))
|
||||
write (777,rec=1) C
|
||||
close(777)
|
||||
call IO_write_jobBinaryFile(777,'F_aim',size(F_aim))
|
||||
write (777,rec=1) F_aim
|
||||
close(777)
|
||||
|
|
|
|||
|
|
@ -123,16 +123,6 @@ program DAMASK_spectral_Driver
|
|||
write(6,'(a)') ' <<<+- DAMASK_spectral_Driver init -+>>>'
|
||||
write(6,'(a)') ' $Id$'
|
||||
#include "compilation_info.f90"
|
||||
write(6,'(a,a)') ' Working Directory: ',trim(getSolverWorkingDirectoryName())
|
||||
write(6,'(a,a)') ' Solver Job Name: ',trim(getSolverJobName())
|
||||
write(6,'(a)') ''
|
||||
write(6,'(a,a)') ' geometry file: ',trim(geometryFile)
|
||||
write(6,'(a)') ''
|
||||
write(6,'(a,3(i12 ))') ' resolution a b c:', res
|
||||
write(6,'(a,3(f12.5))') ' dimension x y z:', geomdim
|
||||
write(6,'(a,i5)') ' homogenization: ', homog
|
||||
write(6,'(a,a)') '',''
|
||||
write(6,'(a,a)') ' Loadcase file: ',trim(loadCaseFile)
|
||||
write(6,'(a)') ''
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! reading basic information from load case file and allocate data structure containing load cases
|
||||
|
|
|
|||
|
|
@ -21,6 +21,19 @@ module DAMASK_spectral_SolverAL
|
|||
solutionState
|
||||
|
||||
implicit none
|
||||
#ifdef PETSC
|
||||
#include <finclude/petscsys.h>
|
||||
#include <finclude/petscvec.h>
|
||||
#include <finclude/petscdmda.h>
|
||||
#include <finclude/petscis.h>
|
||||
#include <finclude/petscmat.h>
|
||||
#include <finclude/petscksp.h>
|
||||
#include <finclude/petscpc.h>
|
||||
#include <finclude/petscsnes.h>
|
||||
#include <finclude/petscvec.h90>
|
||||
#include <finclude/petscsnes.h90>
|
||||
#endif
|
||||
|
||||
character (len=*), parameter, public :: &
|
||||
DAMASK_spectral_SolverAL_label = 'AL'
|
||||
|
||||
|
|
@ -37,7 +50,7 @@ module DAMASK_spectral_SolverAL
|
|||
! PETSc data
|
||||
SNES, private :: snes
|
||||
DM, private :: da
|
||||
Vec, private :: x,r
|
||||
Vec, private :: x,residual
|
||||
PetscMPIInt, private :: rank
|
||||
integer(pInt), private :: iter
|
||||
PetscInt, private :: xs,xm,gxs,gxm
|
||||
|
|
@ -95,18 +108,20 @@ module DAMASK_spectral_SolverAL
|
|||
use mesh, only: &
|
||||
res, &
|
||||
geomdim
|
||||
|
||||
use math, only: &
|
||||
math_invSym3333
|
||||
|
||||
implicit none
|
||||
integer(pInt) :: i,j,k
|
||||
real(pReal), dimension(3,3) :: temp33_Real
|
||||
PetscErrorCode ierr_psc
|
||||
|
||||
PetscObject dummy
|
||||
|
||||
call Utilities_init()
|
||||
|
||||
write(6,'(a)') ''
|
||||
write(6,'(a)') ' <<<+- DAMASK_spectral_solverAL init -+>>>'
|
||||
write(6,'(a)') ' $Id: DAMASK_spectral_SolverAL.f90 1654 2012-08-03 09:25:48Z MPIE\m.diehl $'
|
||||
#include "compilation_info.f90"
|
||||
#include "compilation_info.f90"
|
||||
write(6,'(a)') ''
|
||||
|
||||
allocate (F ( res(1), res(2),res(3),3,3), source = 0.0_pReal)
|
||||
|
|
@ -157,7 +172,7 @@ module DAMASK_spectral_SolverAL
|
|||
coordinates = 0.0 ! change it later!!!
|
||||
endif
|
||||
|
||||
call constitutiveResponse(coordinates,F,F_lastInc,temperature,0.0_pReal,&
|
||||
call Utilities_constitutiveResponse(coordinates,F,F_lastInc,temperature,0.0_pReal,&
|
||||
P,C,P_av,.false.,math_I3)
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
|
|
@ -187,11 +202,11 @@ module DAMASK_spectral_SolverAL
|
|||
DMDA_STENCIL_BOX,res(1),res(2),res(3),PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE, &
|
||||
18,1,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,da,ierr_psc)
|
||||
call DMCreateGlobalVector(da,x,ierr_psc)
|
||||
call VecDuplicate(x,r,ierr_psc)
|
||||
call DMDASetLocalFunction(da,AL_FormRHS,ierr_psc)
|
||||
call VecDuplicate(x,residual,ierr_psc)
|
||||
!call DMDASetLocalFunction(da,AL_FormRHS,ierr_psc)
|
||||
|
||||
call SNESSetDM(snes,da,ierr_psc)
|
||||
call SNESSetFunction(snes,r,SNESDMDAComputeFunction,da,ierr_psc)
|
||||
call SNESSetFunction(snes,residual,AL_FormRHS,dummy,ierr_psc)
|
||||
call SNESSetConvergenceTest(snes,AL_converged,dummy,PETSC_NULL_FUNCTION,ierr_psc)
|
||||
call PetscOptionsInsertString(PetSc_options,ierr_psc)
|
||||
call SNESSetFromOptions(snes,ierr_psc)
|
||||
|
|
@ -239,17 +254,11 @@ module DAMASK_spectral_SolverAL
|
|||
real(pReal), intent(in) :: timeinc, timeinc_old, temperature_bc, guessmode
|
||||
type(boundaryCondition), intent(in) :: P_BC,F_BC
|
||||
real(pReal), dimension(3,3), intent(in) :: rotation_BC
|
||||
|
||||
|
||||
|
||||
SNESConvergedReason reason
|
||||
real(pReal), dimension(3,3) :: deltaF_aim, &
|
||||
F_aim_lab, &
|
||||
F_aim_lab_lastIter
|
||||
F_aim_lab
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! loop variables, convergence etc.
|
||||
real(pReal) :: err_div, err_stress
|
||||
integer(pInt) :: iter, row, column, i, j, k
|
||||
real(pReal) :: defgradDet, defgradDetMax, defgradDetMin
|
||||
real(pReal), dimension(3,3) :: temp33_Real
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
|
|
@ -268,7 +277,7 @@ module DAMASK_spectral_SolverAL
|
|||
write (777,rec=1) C
|
||||
close(777)
|
||||
endif
|
||||
|
||||
AL_solution%converged =.false.
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! winding forward of deformation aim in loadcase system
|
||||
if (F_BC%myType=='l') then ! calculate deltaF_aim from given L and current F
|
||||
|
|
@ -306,7 +315,10 @@ module DAMASK_spectral_SolverAL
|
|||
call AL_InitialGuess(xx_psc)
|
||||
call VecRestoreArrayF90(x,xx_psc,ierr_psc)
|
||||
call SNESSolve(snes,PETSC_NULL_OBJECT,x,ierr_psc)
|
||||
|
||||
call SNESGetConvergedReason(snes,reason,ierr_psc)
|
||||
print*, 'reason', reason
|
||||
if (reason > 0 ) AL_solution%converged = .true.
|
||||
|
||||
end function AL_solution
|
||||
|
||||
! -------------------------------------------------------------------
|
||||
|
|
@ -314,7 +326,7 @@ module DAMASK_spectral_SolverAL
|
|||
subroutine AL_InitialGuess(xx_psc)
|
||||
|
||||
implicit none
|
||||
#include <finclude/petsc.h>
|
||||
|
||||
|
||||
! Input/output variables:
|
||||
|
||||
|
|
@ -347,6 +359,55 @@ module DAMASK_spectral_SolverAL
|
|||
|
||||
return
|
||||
end subroutine AL_InitialGuess
|
||||
|
||||
subroutine AL_FormRHS(snes_local,X_local,F_local,dummy,ierr_psc)
|
||||
|
||||
! Input/output variables:
|
||||
SNES snes_local
|
||||
Vec X_local,F_local
|
||||
PetscErrorCode ierr_psc
|
||||
PetscObject dummy
|
||||
DM da_local
|
||||
|
||||
! Declarations for use with local arrays:
|
||||
PetscScalar,pointer :: lx_v(:),lf_v(:)
|
||||
Vec localX
|
||||
|
||||
! Scatter ghost points to local vector, using the 2-step process
|
||||
! DMGlobalToLocalBegin(), DMGlobalToLocalEnd().
|
||||
! By placing code between these two statements, computations can
|
||||
! be done while messages are in transition.
|
||||
call SNESGetDM(snes_local,da_local,ierr_psc)
|
||||
call DMGetLocalVector(da_local,localX,ierr_psc)
|
||||
call DMGlobalToLocalBegin(da_local,X_local,INSERT_VALUES, &
|
||||
& localX,ierr_psc)
|
||||
call DMGlobalToLocalEnd(da_local,X_local,INSERT_VALUES,localX,ierr_psc)
|
||||
|
||||
! Get a pointer to vector data.
|
||||
! - For default PETSc vectors, VecGetArray90() returns a pointer to
|
||||
! the data array. Otherwise, the routine is implementation dependent.
|
||||
! - You MUST call VecRestoreArrayF90() when you no longer need access to
|
||||
! the array.
|
||||
! - Note that the interface to VecGetArrayF90() differs from VecGetArray(),
|
||||
! and is useable from Fortran-90 Only.
|
||||
|
||||
call VecGetArrayF90(localX,lx_v,ierr_psc)
|
||||
call VecGetArrayF90(F_local,lf_v,ierr_psc)
|
||||
|
||||
! Compute function over the locally owned part of the grid
|
||||
call AL_FormRHS_local(lx_v,lf_v,dummy,ierr_psc)
|
||||
|
||||
! Restore vectors
|
||||
call VecRestoreArrayF90(localX,lx_v,ierr_psc)
|
||||
call VecRestoreArrayF90(F_local,lf_v,ierr_psc)
|
||||
|
||||
! Insert values into global vector
|
||||
|
||||
call DMRestoreLocalVector(da_local,localX,ierr_psc)
|
||||
|
||||
return
|
||||
|
||||
end subroutine AL_FormRHS
|
||||
|
||||
! ---------------------------------------------------------------------
|
||||
!
|
||||
|
|
@ -360,7 +421,7 @@ module DAMASK_spectral_SolverAL
|
|||
! Notes:
|
||||
! This routine uses standard Fortran-style computations over a 3-dim array.
|
||||
!
|
||||
subroutine AL_FormRHS(in,x_scal,f_scal,dummy,ierr_psc)
|
||||
subroutine AL_FormRHS_local(x_scal,f_scal,dummy,ierr_psc)
|
||||
|
||||
use numerics, only: &
|
||||
itmax, &
|
||||
|
|
@ -370,7 +431,8 @@ module DAMASK_spectral_SolverAL
|
|||
math_transpose33, &
|
||||
math_mul3333xx33
|
||||
use mesh, only: &
|
||||
res
|
||||
res, &
|
||||
wgt
|
||||
use DAMASK_spectral_Utilities, only: &
|
||||
field_real, &
|
||||
Utilities_forwardFFT, &
|
||||
|
|
@ -379,13 +441,10 @@ module DAMASK_spectral_SolverAL
|
|||
Utilities_constitutiveResponse
|
||||
|
||||
implicit none
|
||||
#include <finclude/petsc.h>
|
||||
|
||||
integer(pInt) :: i,j,k
|
||||
Input/output variables:
|
||||
DMDALocalInfo in(DMDA_LOCAL_INFO_SIZE)
|
||||
PetscScalar x_scal(0:17,XG_RANGE,YG_RANGE,ZG_RANGE)
|
||||
PetscScalar f_scal(0:17,X_RANGE,Y_RANGE,Z_RANGE)
|
||||
PetscScalar x_scal(0:17,gxs:gxs+gxm,gys:gys+gym,gzs:gzs+gzm)
|
||||
PetscScalar f_scal(0:17,xs:xs+xm,ys:ys+ym,zs:zs+zm)
|
||||
real(pReal), dimension (3,3) :: temp33_real
|
||||
PetscObject dummy
|
||||
PetscErrorCode ierr_psc
|
||||
|
|
@ -399,7 +458,7 @@ module DAMASK_spectral_SolverAL
|
|||
write(6,'(3(a,i6.6))') ' @ Iter. ',itmin,' < ',iter,' < ',itmax
|
||||
write(6,'(a,/,3(3(f12.7,1x)/))',advance='no') 'deformation gradient aim =',&
|
||||
math_transpose33(F_aim)
|
||||
|
||||
|
||||
do k=gzs,gzs+gzm; do j=gys,gys+gym; do i=gxs,gxs+gxm
|
||||
F(i,j,k,1,1) = x_scal(0,i,j,k)
|
||||
F(i,j,k,1,2) = x_scal(1,i,j,k)
|
||||
|
|
@ -420,10 +479,10 @@ module DAMASK_spectral_SolverAL
|
|||
F_lambda(i,j,k,3,2) = x_scal(16,i,j,k)
|
||||
F_lambda(i,j,k,3,3) = x_scal(17,i,j,k)
|
||||
enddo; enddo; enddo
|
||||
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! evaluate constitutive response
|
||||
call constitutiveResponse(coordinates,F,F_lastInc,temperature,params%timeinc,&
|
||||
call Utilities_constitutiveResponse(coordinates,F_lastInc,F,temperature,params%timeinc,&
|
||||
P,C,P_av,ForwardData,params%rotation_BC)
|
||||
ForwardData = .False.
|
||||
|
||||
|
|
@ -432,8 +491,6 @@ module DAMASK_spectral_SolverAL
|
|||
F_aim = F_aim - math_mul3333xx33(S, ((P_av - params%P_BC))) !S = 0.0 for no bc
|
||||
err_stress = maxval(mask_stress * (P_av - params%P_BC)) ! mask = 0.0 for no bc
|
||||
|
||||
|
||||
F_aim_lab = math_rotate_backward33(F_aim,params%rotation_BC) ! boundary conditions from load frame into lab (Fourier) frame
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! doing Fourier transform
|
||||
|
|
@ -441,27 +498,28 @@ module DAMASK_spectral_SolverAL
|
|||
do k = 1_pInt, res(3); do j = 1_pInt, res(2); do i = 1_pInt, res(1)
|
||||
field_real(i,j,k,1:3,1:3) = math_mul3333xx33(C_scale,F_lambda(i,j,k,1:3,1:3)-F(i,j,k,1:3,1:3))
|
||||
enddo; enddo; enddo
|
||||
|
||||
|
||||
|
||||
call Utilities_forwardFFT()
|
||||
call Utilities_fourierConvolution(F_aim_lab)
|
||||
call Utilities_fourierConvolution(math_rotate_backward33(F_aim,params%rotation_BC))
|
||||
call Utilities_backwardFFT()
|
||||
|
||||
err_f = 0.0_pReal
|
||||
err_p = 0.0_pReal
|
||||
|
||||
do k = 1_pInt, res(3); do j = 1_pInt, res(2); do i = 1_pInt, res(1)
|
||||
temp33_real = field_real(i,j,k,1:3,1:3) - F(i,j,k,1:3,1:3)
|
||||
err_f = err_f + sum(temp33_real*temp33_real)
|
||||
|
||||
temp33_real = F_lambda(i,j,k,1:3,1:3) - &
|
||||
math_mul3333xx33(S_scale,P(i,j,k,1:3,1:3)) + math_I3
|
||||
(math_mul3333xx33(S_scale,P(i,j,k,1:3,1:3)) + math_I3)
|
||||
err_p = err_p + sum(temp33_real*temp33_real)
|
||||
enddo; enddo; enddo
|
||||
|
||||
err_f = wgt*sqrt(err_f)/sum((F_aim-math_I3)*(F_aim-math_I3)))
|
||||
err_p = wgt*sqrt(err_p)/sum((F_aim-math_I3)*(F_aim-math_I3)))
|
||||
|
||||
do k=zs,ze; do j=ys,ye; do i=xs,xe
|
||||
enddo; enddo; enddo
|
||||
err_f = wgt*sqrt(err_f)
|
||||
err_p = wgt*sqrt(err_p)
|
||||
|
||||
|
||||
do k=gzs,gzs+gzm; do j=gys,gys+gym; do i=gxs,gxs+gxm
|
||||
temp33_real = math_mul3333xx33(S_scale,P(i,j,k,1:3,1:3)) + math_I3 - F_lambda(i,j,k,1:3,1:3) &
|
||||
+ F(i,j,k,1:3,1:3) - field_real(i,j,k,1:3,1:3)
|
||||
f_scal(0,i,j,k) = temp33_real(1,1)
|
||||
|
|
@ -485,12 +543,12 @@ module DAMASK_spectral_SolverAL
|
|||
enddo; enddo; enddo
|
||||
|
||||
return
|
||||
end subroutine AL_FormRHS
|
||||
end subroutine AL_FormRHS_local
|
||||
|
||||
! ---------------------------------------------------------------------
|
||||
! User defined convergence check
|
||||
!
|
||||
subroutine AL_converged(snes,it,xnorm,snorm,fnorm,reason,dummy,ierr_psc)
|
||||
subroutine AL_converged(snes_local,it,xnorm,snorm,fnorm,reason,dummy,ierr_psc)
|
||||
|
||||
use numerics, only: &
|
||||
itmax, &
|
||||
|
|
@ -501,10 +559,9 @@ module DAMASK_spectral_SolverAL
|
|||
err_stress_tolabs
|
||||
|
||||
implicit none
|
||||
#include <finclude/petsc.h>
|
||||
|
||||
! Input/output variables:
|
||||
SNES snes
|
||||
SNES snes_local
|
||||
PetscInt it
|
||||
PetscReal xnorm, snorm, fnorm
|
||||
SNESConvergedReason reason
|
||||
|
|
@ -512,13 +569,15 @@ module DAMASK_spectral_SolverAL
|
|||
PetscErrorCode ierr_psc
|
||||
logical :: Converged
|
||||
|
||||
Converged = (iter < itmax) .and. (iter > itmin) .and. &
|
||||
Converged = (iter > itmin .and. &
|
||||
all([ err_f/sqrt(sum((F_aim-math_I3)*(F_aim-math_I3)))/err_f_tol, &
|
||||
err_p/sqrt(sum((F_aim-math_I3)*(F_aim-math_I3)))/err_p_tol, &
|
||||
err_stress/min(maxval(abs(P_av))*err_stress_tolrel,err_stress_tolabs)] < 1.0_pReal)
|
||||
err_stress/min(maxval(abs(P_av))*err_stress_tolrel,err_stress_tolabs)] < 1.0_pReal))
|
||||
|
||||
if (Converged) then
|
||||
reason = 1
|
||||
elseif (iter > itmax) then
|
||||
reason = -1
|
||||
else
|
||||
reason = 0
|
||||
endif
|
||||
|
|
@ -527,14 +586,17 @@ module DAMASK_spectral_SolverAL
|
|||
write(6,'(a,es14.7)') 'error F = ', err_f/sqrt(sum((F_aim-math_I3)*(F_aim-math_I3)))/err_f_tol
|
||||
write(6,'(a,es14.7)') 'error P = ', err_p/sqrt(sum((F_aim-math_I3)*(F_aim-math_I3)))/err_p_tol
|
||||
return
|
||||
|
||||
end subroutine AL_converged
|
||||
|
||||
subroutine AL_destroy()
|
||||
|
||||
use DAMASK_spectral_Utilities, only: &
|
||||
Utilities_destroy
|
||||
implicit none
|
||||
PetscErrorCode ierr_psc
|
||||
|
||||
call VecDestroy(x,ierr_psc)
|
||||
call VecDestroy(r,ierr_psc)
|
||||
call VecDestroy(residual,ierr_psc)
|
||||
call SNESDestroy(snes,ierr_psc)
|
||||
call DMDestroy(da,ierr_psc)
|
||||
call PetscFinalize(ierr_psc)
|
||||
|
|
|
|||
|
|
@ -178,15 +178,19 @@ subroutine DAMASK_interface_init(loadCaseParameterIn,geometryParameterIn)
|
|||
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()
|
||||
write(6,'(a,a)') ' Command line call: ', trim(commandLine)
|
||||
write(6,'(a,a)') ' Geometry parameter: ', trim(geometryParameter)
|
||||
write(6,'(a,a)') ' Loadcase parameter: ', trim(loadcaseParameter)
|
||||
write(6,'(a,a)') ' Host name: ', trim(hostName)
|
||||
write(6,'(a,a)') ' User name: ', trim(userName)
|
||||
write(6,'(a,a)') ' Path separator: ', getPathSep()
|
||||
write(6,'(a,a)') ' Command line call: ', trim(commandLine)
|
||||
write(6,'(a,a)') ' Geometry parameter: ', trim(geometryParameter)
|
||||
write(6,'(a,a)') ' Loadcase parameter: ', trim(loadcaseParameter)
|
||||
write(6,'(a,a)') ' Geometry file: ', trim(geometryFile)
|
||||
write(6,'(a,a)') ' Loadcase file: ', trim(loadCaseFile)
|
||||
write(6,'(a,a)') ' Working Directory: ', trim(getSolverWorkingDirectoryName())
|
||||
write(6,'(a,a)') ' Solver Job Name: ', trim(getSolverJobName())
|
||||
if (SpectralRestartInc > 1_pInt) write(6,'(a,i6.6)') &
|
||||
' Restart at increment: ', spectralRestartInc
|
||||
write(6,'(a,l1)') ' Append to result file: ', appendToOutFile
|
||||
' Restart at increment: ', spectralRestartInc
|
||||
write(6,'(a,l1,/)') ' Append to result file: ', appendToOutFile
|
||||
|
||||
end subroutine DAMASK_interface_init
|
||||
|
||||
|
|
|
|||
|
|
@ -100,9 +100,14 @@ subroutine FE_init
|
|||
character(len=1024) :: line
|
||||
integer(pInt), dimension(1_pInt+2_pInt*maxNchunks) :: positions
|
||||
#endif
|
||||
!$OMP CRITICAL (write2out)
|
||||
write(6,*)
|
||||
write(6,*) '<<<+- FEsolving init -+>>>'
|
||||
write(6,*) '$Id$'
|
||||
#include "compilation_info.f90"
|
||||
!$OMP END CRITICAL (write2out)
|
||||
|
||||
modelName = getSolverJobName()
|
||||
|
||||
#ifdef Spectral
|
||||
restartInc = spectralRestartInc
|
||||
if(restartInc <= 0_pInt) then
|
||||
|
|
@ -169,13 +174,11 @@ subroutine FE_init
|
|||
#endif
|
||||
200 close(fileunit)
|
||||
endif
|
||||
|
||||
! the following array are allocated by mesh.f90 and need to be deallocated in case of regridding
|
||||
if (allocated(calcMode)) deallocate(calcMode)
|
||||
if (allocated(FEsolving_execIP)) deallocate(FEsolving_execIP)
|
||||
#endif
|
||||
!$OMP CRITICAL (write2out)
|
||||
write(6,*)
|
||||
write(6,*) '<<<+- FEsolving init -+>>>'
|
||||
write(6,*) '$Id$'
|
||||
#include "compilation_info.f90"
|
||||
if (iand(debug_level(debug_FEsolving),debug_levelBasic) /= 0_pInt) then
|
||||
write(6,*) 'restart writing: ', restartWrite
|
||||
write(6,*) 'restart reading: ', restartRead
|
||||
|
|
|
|||
|
|
@ -39,10 +39,10 @@ SHELL = /bin/sh
|
|||
########################################################################################
|
||||
|
||||
#auto values will be set by setup_code.py
|
||||
FFTWROOT := auto
|
||||
IMKLROOT := auto
|
||||
ACMLROOT := auto
|
||||
LAPACKROOT := auto
|
||||
FFTWROOT := /usr/local
|
||||
IMKLROOT :=
|
||||
ACMLROOT :=
|
||||
LAPACKROOT := /usr
|
||||
|
||||
F90 ?= ifort
|
||||
COMPILERNAME ?= $(F90)
|
||||
|
|
@ -96,11 +96,13 @@ endif
|
|||
LIBRARIES +=-lfftw3
|
||||
LIB_DIRS +=-L$(FFTWROOT)/lib
|
||||
|
||||
#ifdef PETSC_DIR
|
||||
#include ${PETSC_DIR}/conf/variables
|
||||
#INCLUDE_DIRS +=${PETSC_FC_INCLUDES} -DPETSC
|
||||
#LIBRARIES +=${PETSC_WITH_EXTERNAL_LIB}
|
||||
#endif
|
||||
ifeq "$(SOLVER)" "NEW"
|
||||
ifdef PETSC_DIR
|
||||
include ${PETSC_DIR}/conf/variables
|
||||
INCLUDE_DIRS +=${PETSC_FC_INCLUDES} -DPETSC
|
||||
LIBRARIES +=${PETSC_WITH_EXTERNAL_LIB}
|
||||
endif
|
||||
endif
|
||||
|
||||
ifdef IMKLROOT
|
||||
LIB_DIRS +=-L$(IMKLROOT)/lib/intel64
|
||||
|
|
@ -273,10 +275,29 @@ COMPILE_MAXOPTI =$(OPENMP_FLAG_$(F90)) $(COMPILE_OPTIONS_$(F90)) $(STANDARD_CHEC
|
|||
COMPILED_FILES = prec.o DAMASK_spectral_interface.o IO.o numerics.o debug.o math.o \
|
||||
FEsolving.o mesh.o material.o lattice.o \
|
||||
constitutive_dislotwin.o constitutive_j2.o constitutive_phenopowerlaw.o \
|
||||
constitutive_titanmod.o constitutive_nonlocal.o constitutive_none.o constitutive.o \
|
||||
homogenization_RGC.o homogenization_isostrain.o homogenization.o CPFEM.o crystallite.o
|
||||
constitutive_titanmod.o constitutive_nonlocal.o constitutive_none.o constitutive.o crystallite.o \
|
||||
homogenization_RGC.o homogenization_isostrain.o homogenization.o CPFEM.o
|
||||
|
||||
ifneq "$(SOLVER)" "AL"
|
||||
ifeq "$(SOLVER)" "NEW"
|
||||
COMPILED_FILES += DAMASK_spectral_Utilities.o DAMASK_spectral_SolverBasic.o DAMASK_spectral_SolverAL.o
|
||||
|
||||
DAMASK_spectral.exe: DAMASK_spectral_Driver.o
|
||||
$(PREFIX) $(COMPILERNAME) $(OPENMP_FLAG_$(F90)) $(OPTIMIZATION_$(MAXOPTI)_$(F90)) $(STANDARD_CHECK_$(F90)) \
|
||||
-o DAMASK_spectral.exe DAMASK_spectral_Driver.o \
|
||||
$(COMPILED_FILES) $(LIB_DIRS) $(LIBRARIES) $(SUFFIX)
|
||||
|
||||
DAMASK_spectral_Driver.o: DAMASK_spectral_Driver.f90 DAMASK_spectral_SolverBasic.o DAMASK_spectral_SolverAL.o
|
||||
$(PREFIX) $(COMPILERNAME) $(COMPILE_MAXOPTI) -c DAMASK_spectral_Driver.f90 $(SUFFIX)
|
||||
|
||||
DAMASK_spectral_SolverAL.o: DAMASK_spectral_SolverAL.f90\
|
||||
DAMASK_spectral_Utilities.o
|
||||
|
||||
DAMASK_spectral_SolverBasic.o: DAMASK_spectral_SolverBasic.f90\
|
||||
DAMASK_spectral_Utilities.o
|
||||
|
||||
DAMASK_spectral_Utilities.o: DAMASK_spectral_Utilities.f90\
|
||||
CPFEM.o
|
||||
else
|
||||
DAMASK_spectral.exe: DAMASK_spectral.o
|
||||
$(PREFIX) $(COMPILERNAME) $(OPENMP_FLAG_$(F90)) $(OPTIMIZATION_$(MAXOPTI)_$(F90)) $(STANDARD_CHECK_$(F90)) \
|
||||
-o DAMASK_spectral.exe DAMASK_spectral.o \
|
||||
|
|
@ -284,15 +305,7 @@ DAMASK_spectral.exe: DAMASK_spectral.o
|
|||
|
||||
DAMASK_spectral.o: DAMASK_spectral.f90 CPFEM.o
|
||||
$(PREFIX) $(COMPILERNAME) $(COMPILE_MAXOPTI) -c DAMASK_spectral.f90 $(SUFFIX)
|
||||
else
|
||||
DAMASK_spectral_AL.exe: DAMASK_spectral_AL.o
|
||||
$(PREFIX) $(COMPILERNAME) $(OPENMP_FLAG_$(F90)) $(OPTIMIZATION_$(MAXOPTI)_$(F90)) $(STANDARD_CHECK_$(F90)) \
|
||||
-o DAMASK_spectral_AL.exe DAMASK_spectral_AL.o \
|
||||
$(COMPILED_FILES) $(LIB_DIRS) $(LIBRARIES) $(SUFFIX)
|
||||
|
||||
DAMASK_spectral_AL.o: DAMASK_spectral_AL.f90 CPFEM.o
|
||||
$(PREFIX) $(COMPILERNAME) $(COMPILE_MAXOPTI) -c DAMASK_spectral_AL.f90 $(SUFFIX)
|
||||
endif
|
||||
endif
|
||||
CPFEM.o: CPFEM.f90\
|
||||
homogenization.o
|
||||
|
||||
|
|
|
|||
143
code/mesh.f90
143
code/mesh.f90
|
|
@ -29,7 +29,7 @@
|
|||
|
||||
module mesh
|
||||
use prec, only: pReal, pInt
|
||||
|
||||
|
||||
implicit none
|
||||
private
|
||||
|
||||
|
|
@ -94,7 +94,14 @@ module mesh
|
|||
mesh_subNodeCoord !< coordinates of subnodes per element
|
||||
|
||||
logical, private :: noPart !< for cases where the ABAQUS input file does not use part/assembly information
|
||||
|
||||
|
||||
#ifdef Spectral
|
||||
real(pReal), dimension(3), public :: geomdim, virt_dim ! physical dimension of volume element per direction
|
||||
integer(pInt), dimension(3), public :: res
|
||||
real(pReal), public :: wgt
|
||||
integer(pInt), public :: res1_red, homog
|
||||
integer(pInt), private :: i
|
||||
#endif
|
||||
|
||||
! These definitions should actually reside in the FE-solver specific part (different for MARC/ABAQUS)
|
||||
! Hence, I suggest to prefix with "FE_"
|
||||
|
|
@ -272,14 +279,14 @@ module mesh
|
|||
mesh_build_ipVolumes, &
|
||||
mesh_build_ipCoordinates
|
||||
#ifdef Spectral
|
||||
public :: mesh_spectral_getResolution, &
|
||||
mesh_spectral_getDimension, &
|
||||
mesh_spectral_getHomogenization, &
|
||||
mesh_regrid
|
||||
public :: mesh_regrid
|
||||
#endif
|
||||
|
||||
private :: &
|
||||
#ifdef Spectral
|
||||
mesh_spectral_getResolution, &
|
||||
mesh_spectral_getDimension, &
|
||||
mesh_spectral_getHomogenization, &
|
||||
mesh_spectral_count_nodesAndElements, &
|
||||
mesh_spectral_count_cpElements, &
|
||||
mesh_spectral_map_elements, &
|
||||
|
|
@ -340,6 +347,8 @@ subroutine mesh_init(ip,element)
|
|||
#endif
|
||||
#ifdef Spectral
|
||||
IO_open_file
|
||||
use numerics, only: &
|
||||
divergence_correction
|
||||
#else
|
||||
IO_open_InputFile
|
||||
#endif
|
||||
|
|
@ -375,8 +384,6 @@ subroutine mesh_init(ip,element)
|
|||
if (allocated(mesh_ipNeighborhood)) deallocate(mesh_ipNeighborhood)
|
||||
if (allocated(mesh_ipVolume)) deallocate(mesh_ipVolume)
|
||||
if (allocated(mesh_nodeTwins)) deallocate(mesh_nodeTwins)
|
||||
if (allocated(calcMode)) deallocate(calcMode)
|
||||
if (allocated(FEsolving_execIP)) deallocate(FEsolving_execIP)
|
||||
if (allocated(FE_nodesAtIP)) deallocate(FE_nodesAtIP)
|
||||
if (allocated(FE_ipNeighbor))deallocate(FE_ipNeighbor)
|
||||
if (allocated(FE_subNodeParent)) deallocate(FE_subNodeParent)
|
||||
|
|
@ -384,12 +391,27 @@ subroutine mesh_init(ip,element)
|
|||
call mesh_build_FEdata ! get properties of the different types of elements
|
||||
#ifdef Spectral
|
||||
call IO_open_file(fileUnit,geometryFile) ! parse info from geometry file...
|
||||
call mesh_spectral_count_nodesAndElements(fileUnit)
|
||||
res = mesh_spectral_getResolution(fileUnit)
|
||||
res1_red = res(1)/2_pInt + 1_pInt
|
||||
wgt = 1.0/real(res(1)*res(2)*res(3),pReal)
|
||||
geomdim = mesh_spectral_getDimension(fileUnit)
|
||||
homog = mesh_spectral_getHomogenization(fileUnit)
|
||||
if (divergence_correction) then
|
||||
do i = 1_pInt, 3_pInt
|
||||
if (i /= minloc(geomdim,1) .and. i /= maxloc(geomdim,1)) virt_dim = geomdim/geomdim(i)
|
||||
enddo
|
||||
else
|
||||
virt_dim = geomdim
|
||||
endif
|
||||
write(6,'(a,3(i12 ))') ' resolution a b c:', res
|
||||
write(6,'(a,3(f12.5))') ' dimension x y z:', geomdim
|
||||
write(6,'(a,i5,/)') ' homogenization: ', homog
|
||||
call mesh_spectral_count_nodesAndElements
|
||||
call mesh_spectral_count_cpElements
|
||||
call mesh_spectral_map_elements
|
||||
call mesh_spectral_map_nodes
|
||||
call mesh_spectral_count_cpSizes
|
||||
call mesh_spectral_build_nodes(fileUnit)
|
||||
call mesh_spectral_build_nodes
|
||||
call mesh_spectral_build_elements(fileUnit)
|
||||
#endif
|
||||
#ifdef Marc
|
||||
|
|
@ -864,10 +886,10 @@ function mesh_regrid(adaptive,resNewInput,minRes)
|
|||
integer(pInt), dimension(3), optional, intent(in) :: minRes
|
||||
integer(pInt), dimension(3) :: mesh_regrid, ratio
|
||||
integer(pInt), dimension(3,2) :: possibleResNew
|
||||
integer(pInt):: maxsize, i, j, k, ielem, Npoints, NpointsNew, spatialDim
|
||||
integer(pInt), dimension(3) :: res, resNew
|
||||
integer(pInt):: maxsize, i, j, k, ielem, NpointsNew, spatialDim
|
||||
integer(pInt), dimension(3) :: resNew
|
||||
integer(pInt), dimension(:), allocatable :: indices
|
||||
real(pReal), dimension(3) :: geomdim,geomdimNew
|
||||
real(pReal), dimension(3) :: geomdimNew
|
||||
real(pReal), dimension(3,3) :: Favg, Favg_LastInc
|
||||
real(pReal), dimension(:,:), allocatable :: &
|
||||
coordinatesNew, &
|
||||
|
|
@ -878,10 +900,10 @@ function mesh_regrid(adaptive,resNewInput,minRes)
|
|||
|
||||
real(pReal), dimension(:,:,:,:,:), allocatable :: &
|
||||
F, FNew, &
|
||||
F_lastInc, F_lastIncNew, &
|
||||
Fp, FpNew, &
|
||||
Lp, LpNew, &
|
||||
dcsdE, dcsdENew
|
||||
dcsdE, dcsdENew, &
|
||||
F_lastIncNew
|
||||
|
||||
real(pReal), dimension (:,:,:,:,:,:,:), allocatable :: &
|
||||
dPdF, dPdFNew
|
||||
|
|
@ -892,10 +914,6 @@ function mesh_regrid(adaptive,resNewInput,minRes)
|
|||
stateConst
|
||||
integer(pInt), dimension(:,:), allocatable :: &
|
||||
sizeStateConst, sizeStateHomog
|
||||
|
||||
res = mesh_spectral_getResolution()
|
||||
geomdim = mesh_spectral_getDimension()
|
||||
Npoints = res(1)*res(2)*res(3)
|
||||
|
||||
if (adaptive) then
|
||||
if (present(resNewInput)) then
|
||||
|
|
@ -911,16 +929,17 @@ function mesh_regrid(adaptive,resNewInput,minRes)
|
|||
read (777,rec=1) F
|
||||
close (777)
|
||||
|
||||
! ----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)) / real(Npoints,pReal)
|
||||
enddo; enddo
|
||||
! ----read in average deformation--------
|
||||
call IO_read_jobBinaryFile(777,'F_aim',trim(getSolverJobName()),size(Favg))
|
||||
read (777,rec=1) Favg
|
||||
close (777)
|
||||
|
||||
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))
|
||||
allocate(coordinatesLinear(3,mesh_NcpElems))
|
||||
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
|
||||
|
|
@ -954,7 +973,6 @@ function mesh_regrid(adaptive,resNewInput,minRes)
|
|||
ratio = floor(real(resNewInput,pReal) * (geomdimNew/geomdim), pInt)
|
||||
|
||||
possibleResNew = 1_pInt
|
||||
|
||||
do i = 1_pInt, spatialDim
|
||||
if (mod(ratio(i),2) == 0_pInt) then
|
||||
possibleResNew(i,1:2) = [ratio(i),ratio(i) + 2_pInt]
|
||||
|
|
@ -1016,7 +1034,7 @@ function mesh_regrid(adaptive,resNewInput,minRes)
|
|||
|
||||
!----- Nearest neighbour search ------------------------------------
|
||||
allocate(indices(NpointsNew))
|
||||
call math_nearestNeighborSearch(spatialDim, Favg, geomdim, NpointsNew, Npoints, &
|
||||
call math_nearestNeighborSearch(spatialDim, Favg, geomdim, NpointsNew, mesh_NcpElems, &
|
||||
coordinatesNew, coordinatesLinear, indices)
|
||||
deallocate(coordinatesNew)
|
||||
|
||||
|
|
@ -1045,7 +1063,7 @@ function mesh_regrid(adaptive,resNewInput,minRes)
|
|||
|
||||
call IO_write_jobFile(777,'idx2') ! make it a general open-write file
|
||||
write(777, '(A)') '1 header'
|
||||
write(777, '(A)') 'Numbered indices as per the large set'
|
||||
write(777, '(A)') 'Numbered indices as per the small set'
|
||||
do i = 1_pInt, NpointsNew
|
||||
write(777,trim(formatString),advance='no') indices(i), ' '
|
||||
if(mod(i,resNew(1)) == 0_pInt) write(777,'(A)') ''
|
||||
|
|
@ -1053,7 +1071,7 @@ function mesh_regrid(adaptive,resNewInput,minRes)
|
|||
close(777)
|
||||
|
||||
!------Adjusting the point-to-grain association---------------------
|
||||
write(N_Digits, '(I16.16)') 1_pInt+int(log10(real(maxval(mesh_element(4,1:Npoints)),pReal)),pInt)
|
||||
write(N_Digits, '(I16.16)') 1_pInt+int(log10(real(maxval(mesh_element(4,1:mesh_NcpElems)),pReal)),pInt)
|
||||
N_Digits = adjustl(N_Digits)
|
||||
formatString = '(I'//trim(N_Digits)//'.'//trim(N_Digits)//',a)'
|
||||
|
||||
|
|
@ -1069,25 +1087,24 @@ function mesh_regrid(adaptive,resNewInput,minRes)
|
|||
enddo
|
||||
close(777)
|
||||
|
||||
!---------------------------------------------------------
|
||||
allocate(F_lastInc(res(1),res(2),res(3),3,3))
|
||||
allocate(F_lastIncNew(resNew(1),resNew(2),resNew(3),3,3))
|
||||
call IO_read_jobBinaryFile(777,'convergedSpectralDefgrad_lastInc',trim(getSolverJobName()),size(F_lastInc))
|
||||
read (777,rec=1) F_lastInc
|
||||
!---set F_lastInc to homogeneous deformation------------------------------------------------------
|
||||
|
||||
call IO_read_jobBinaryFile(777,'F_aim_lastInc',trim(getSolverJobName()),size(Favg_LastInc))
|
||||
read (777,rec=1) Favg_LastInc
|
||||
close (777)
|
||||
do i= 1_pInt,3_pInt; do j = 1_pInt,3_pInt
|
||||
Favg_LastInc(i,j) = sum(F_lastInc(1:res(1),1:res(2),1:res(3),i,j)) / real(Npoints,pReal)
|
||||
enddo; enddo
|
||||
|
||||
allocate(F_lastIncNew(resNew(1),resNew(2),resNew(3),3,3))
|
||||
do k=1_pInt,resNew(3); do j=1_pInt, resNew(2); do i=1_pInt, resNew(1)
|
||||
F_lastIncNew(i,j,k,1:3,1:3) = Favg
|
||||
F_lastIncNew(i,j,k,1:3,1:3) = Favg_LastInc
|
||||
enddo; enddo; enddo
|
||||
|
||||
call IO_write_jobBinaryFile(777,'convergedSpectralDefgrad_lastInc',size(F_lastIncNew))
|
||||
write (777,rec=1) F_lastIncNew
|
||||
close (777)
|
||||
deallocate(F_lastInc)
|
||||
deallocate(F_lastIncNew)
|
||||
|
||||
!-------------------------------------------------------------------
|
||||
allocate(material_phase (1,1, Npoints))
|
||||
allocate(material_phase (1,1, mesh_NcpElems))
|
||||
allocate(material_phaseNew (1,1, NpointsNew))
|
||||
call IO_read_jobBinaryFile(777,'recordedPhase',trim(getSolverJobName()),size(material_phase))
|
||||
read (777,rec=1) material_phase
|
||||
|
|
@ -1103,7 +1120,7 @@ function mesh_regrid(adaptive,resNewInput,minRes)
|
|||
deallocate(material_phaseNew)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
allocate(F (3,3,1,1, Npoints))
|
||||
allocate(F (3,3,1,1, mesh_NcpElems))
|
||||
allocate(FNew (3,3,1,1, NpointsNew))
|
||||
call IO_read_jobBinaryFile(777,'convergedF',trim(getSolverJobName()),size(F))
|
||||
read (777,rec=1) F
|
||||
|
|
@ -1118,7 +1135,7 @@ function mesh_regrid(adaptive,resNewInput,minRes)
|
|||
deallocate(F)
|
||||
deallocate(FNew)
|
||||
!---------------------------------------------------------------------
|
||||
allocate(Fp (3,3,1,1,Npoints))
|
||||
allocate(Fp (3,3,1,1,mesh_NcpElems))
|
||||
allocate(FpNew (3,3,1,1,NpointsNew))
|
||||
call IO_read_jobBinaryFile(777,'convergedFp',trim(getSolverJobName()),size(Fp))
|
||||
read (777,rec=1) Fp
|
||||
|
|
@ -1134,7 +1151,7 @@ function mesh_regrid(adaptive,resNewInput,minRes)
|
|||
deallocate(FpNew)
|
||||
|
||||
!------------------------------------------------------------------------
|
||||
allocate(Lp (3,3,1,1,Npoints))
|
||||
allocate(Lp (3,3,1,1,mesh_NcpElems))
|
||||
allocate(LpNew (3,3,1,1,NpointsNew))
|
||||
call IO_read_jobBinaryFile(777,'convergedLp',trim(getSolverJobName()),size(Lp))
|
||||
read (777,rec=1) Lp
|
||||
|
|
@ -1149,7 +1166,7 @@ function mesh_regrid(adaptive,resNewInput,minRes)
|
|||
deallocate(LpNew)
|
||||
|
||||
!----------------------------------------------------------------------------
|
||||
allocate(dcsdE (6,6,1,1,Npoints))
|
||||
allocate(dcsdE (6,6,1,1,mesh_NcpElems))
|
||||
allocate(dcsdENew (6,6,1,1,NpointsNew))
|
||||
call IO_read_jobBinaryFile(777,'convergeddcsdE',trim(getSolverJobName()),size(dcsdE))
|
||||
read (777,rec=1) dcsdE
|
||||
|
|
@ -1164,7 +1181,7 @@ function mesh_regrid(adaptive,resNewInput,minRes)
|
|||
deallocate(dcsdENew)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
allocate(dPdF (3,3,3,3,1,1,Npoints))
|
||||
allocate(dPdF (3,3,3,3,1,1,mesh_NcpElems))
|
||||
allocate(dPdFNew (3,3,3,3,1,1,NpointsNew))
|
||||
call IO_read_jobBinaryFile(777,'convergeddPdF',trim(getSolverJobName()),size(dPdF))
|
||||
read (777,rec=1) dPdF
|
||||
|
|
@ -1179,7 +1196,7 @@ function mesh_regrid(adaptive,resNewInput,minRes)
|
|||
deallocate(dPdFNew)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
allocate(Tstar (6,1,1,Npoints))
|
||||
allocate(Tstar (6,1,1,mesh_NcpElems))
|
||||
allocate(TstarNew (6,1,1,NpointsNew))
|
||||
call IO_read_jobBinaryFile(777,'convergedTstar',trim(getSolverJobName()),size(Tstar))
|
||||
read (777,rec=1) Tstar
|
||||
|
|
@ -1194,16 +1211,16 @@ function mesh_regrid(adaptive,resNewInput,minRes)
|
|||
deallocate(TstarNew)
|
||||
|
||||
!----------------------------------------------------------------------------
|
||||
allocate(sizeStateConst(1,Npoints))
|
||||
allocate(sizeStateConst(1,mesh_NcpElems))
|
||||
call IO_read_jobBinaryFile(777,'sizeStateConst',trim(getSolverJobName()),size(sizeStateConst))
|
||||
read (777,rec=1) sizeStateConst
|
||||
close (777)
|
||||
maxsize = maxval(sizeStateConst(1,1:Npoints))
|
||||
allocate(StateConst (1,1,Npoints,maxsize))
|
||||
maxsize = maxval(sizeStateConst(1,1:mesh_NcpElems))
|
||||
allocate(StateConst (1,1,mesh_NcpElems,maxsize))
|
||||
|
||||
call IO_read_jobBinaryFile(777,'convergedStateConst',trim(getSolverJobName()))
|
||||
k = 0_pInt
|
||||
do i =1, Npoints
|
||||
do i =1, mesh_NcpElems
|
||||
do j = 1,sizeStateConst(1,i)
|
||||
k = k+1_pInt
|
||||
read(777,rec=k) StateConst(1,1,i,j)
|
||||
|
|
@ -1223,16 +1240,16 @@ function mesh_regrid(adaptive,resNewInput,minRes)
|
|||
deallocate(StateConst)
|
||||
|
||||
!----------------------------------------------------------------------------
|
||||
allocate(sizeStateHomog(1,Npoints))
|
||||
allocate(sizeStateHomog(1,mesh_NcpElems))
|
||||
call IO_read_jobBinaryFile(777,'sizeStateHomog',trim(getSolverJobName()),size(sizeStateHomog))
|
||||
read (777,rec=1) sizeStateHomog
|
||||
close (777)
|
||||
maxsize = maxval(sizeStateHomog(1,1:Npoints))
|
||||
allocate(stateHomog (1,1,Npoints,maxsize))
|
||||
maxsize = maxval(sizeStateHomog(1,1:mesh_NcpElems))
|
||||
allocate(stateHomog (1,1,mesh_NcpElems,maxsize))
|
||||
|
||||
call IO_read_jobBinaryFile(777,'convergedStateHomog',trim(getSolverJobName()))
|
||||
k = 0_pInt
|
||||
do i =1, Npoints
|
||||
do i =1, mesh_NcpElems
|
||||
do j = 1,sizeStateHomog(1,i)
|
||||
k = k+1_pInt
|
||||
read(777,rec=k) stateHomog(1,1,i,j)
|
||||
|
|
@ -1261,13 +1278,9 @@ end function mesh_regrid
|
|||
!> @brief Count overall number of nodes and elements in mesh and stores them in
|
||||
!! 'mesh_Nelems' and 'mesh_Nnodes'
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
subroutine mesh_spectral_count_nodesAndElements(myUnit)
|
||||
subroutine mesh_spectral_count_nodesAndElements()
|
||||
|
||||
implicit none
|
||||
integer(pInt), intent(in) :: myUnit
|
||||
integer(pInt), dimension(3) :: res
|
||||
|
||||
res = mesh_spectral_getResolution(myUnit)
|
||||
mesh_Nelems = res(1)*res(2)*res(3)
|
||||
mesh_Nnodes = (1_pInt + res(1))*(1_pInt + res(2))*(1_pInt + res(3))
|
||||
|
||||
|
|
@ -1344,22 +1357,16 @@ end subroutine mesh_spectral_count_cpSizes
|
|||
!> @brief Store x,y,z coordinates of all nodes in mesh.
|
||||
!! Allocates global arrays 'mesh_node0' and 'mesh_node'
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
subroutine mesh_spectral_build_nodes(myUnit)
|
||||
subroutine mesh_spectral_build_nodes()
|
||||
|
||||
use IO, only: &
|
||||
IO_error
|
||||
|
||||
implicit none
|
||||
integer(pInt), intent(in) :: myUnit
|
||||
integer(pInt) :: n
|
||||
integer(pInt), dimension(3) :: res = 1_pInt
|
||||
real(pReal), dimension(3) :: geomdim = 1.0_pReal
|
||||
|
||||
allocate ( mesh_node0 (3,mesh_Nnodes) ); mesh_node0 = 0.0_pReal
|
||||
allocate ( mesh_node (3,mesh_Nnodes) ); mesh_node = 0.0_pReal
|
||||
|
||||
res = mesh_spectral_getResolution(myUnit)
|
||||
geomdim = mesh_spectral_getDimension(myUnit)
|
||||
|
||||
forall (n = 0_pInt:mesh_Nnodes-1_pInt)
|
||||
mesh_node0(1,n+1_pInt) = &
|
||||
|
|
@ -1399,16 +1406,12 @@ subroutine mesh_spectral_build_elements(myUnit)
|
|||
|
||||
integer(pInt), parameter :: maxNchunks = 7_pInt
|
||||
integer(pInt), dimension (1_pInt+2_pInt*maxNchunks) :: myPos
|
||||
integer(pInt), dimension(3) :: res
|
||||
integer(pInt) :: e, i, homog = 0_pInt, headerLength = 0_pInt, maxIntCount
|
||||
integer(pInt) :: e, i, headerLength = 0_pInt, maxIntCount
|
||||
integer(pInt), dimension(:), allocatable :: microstructures
|
||||
integer(pInt), dimension(1,1) :: dummySet = 0_pInt
|
||||
character(len=65536) :: line,keyword
|
||||
character(len=64), dimension(1) :: dummyName = ''
|
||||
|
||||
res = mesh_spectral_getResolution(myUnit)
|
||||
homog = mesh_spectral_getHomogenization(myUnit)
|
||||
|
||||
call IO_checkAndRewind(myUnit)
|
||||
|
||||
read(myUnit,'(a65536)') line
|
||||
|
|
|
|||
|
|
@ -419,7 +419,7 @@ subroutine numerics_init
|
|||
.not. memory_efficient) call IO_error(error_ID = 847_pInt)
|
||||
#endif
|
||||
if (fixedSeed <= 0_pInt) then
|
||||
write(6,'(a)') ' Random is random!'
|
||||
write(6,'(a,/)') ' Random is random!'
|
||||
endif
|
||||
|
||||
end subroutine numerics_init
|
||||
|
|
|
|||
Loading…
Reference in New Issue