223031012
/
DAMASK_EICMD
mirror of https://github.com/achalhp/DAMASK_EICMD.git
1
0
Fork 0
Code Issues Projects Releases Wiki Activity
DAMASK_EICMD/code/DAMASK_spectral_solverAL.f90

697 lines
31 KiB
Fortran
Raw Normal View History

updated copyright header
2013-03-22 23:05:05 +05:30
! Copyright 2011-13 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/>.
!
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
!--------------------------------------------------------------------------------------------------
added svn properties
2013-08-01 14:43:46 +05:30
! $Id$
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
!--------------------------------------------------------------------------------------------------
!> @author Pratheek Shanthraj, Max-Planck-Institut für Eisenforschung GmbH
!> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH
!> @author Philip Eisenlohr, Max-Planck-Institut für Eisenforschung GmbH
!> @brief AL scheme solver
!--------------------------------------------------------------------------------------------------
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
module DAMASK_spectral_solverAL
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
use prec, only: &
pInt, &
pReal
use math, only: &
math_I3
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
use DAMASK_spectral_utilities, only: &
numerics got some protected statements and is now reading in the new keywords for the solver selection in small letters. spectral solver got cut back facilities + improved output to screen
2012-10-02 20:56:56 +05:30
tSolutionState
next round in modularization
2012-07-25 19:31:39 +05:30
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
implicit none
added missing line continuation in numerics and explicit private statements in the other files
2013-03-28 16:07:00 +05:30
private
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)
2012-08-06 18:13:05 +05:30
#include <finclude/petscsys.h>
#include <finclude/petscdmda.h>
#include <finclude/petscsnes.h>
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
character (len=*), parameter, public :: &
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
DAMASK_spectral_solverAL_label = 'al'
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
!--------------------------------------------------------------------------------------------------
introduced possibility so specify deformation gradient aim at end of load case, rate will be calculated using difference between start of load case and aim. Needed for cyclic loading. Use keyword "f" for this behavior, don't use it as short name for freq any more!
2013-05-13 15:14:23 +05:30
! derived types
made doxygen working for MSC.Marc again, small changes on the todo-statements for doxygen, they don't work on single lines of codes but only on module/variables
2013-06-11 22:05:04 +05:30
type tSolutionParams !< @todo use here the type definition for a full loadcase including mask
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
real(pReal), dimension(3,3) :: P_BC, rotation_BC
real(pReal) :: timeinc
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
real(pReal) :: timeincOld
fixed wrong temperature when using spectral solver
2013-01-09 23:38:08 +05:30
real(pReal) :: temperature
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
real(pReal) :: density
numerics got some protected statements and is now reading in the new keywords for the solver selection in small letters. spectral solver got cut back facilities + improved output to screen
2012-10-02 20:56:56 +05:30
end type tSolutionParams
more modularization….
2012-07-26 19:28:47 +05:30
numerics got some protected statements and is now reading in the new keywords for the solver selection in small letters. spectral solver got cut back facilities + improved output to screen
2012-10-02 20:56:56 +05:30
type(tSolutionParams), private :: params
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
real(pReal), private, dimension(3,3) :: mask_stress = 0.0_pReal
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
!--------------------------------------------------------------------------------------------------
! PETSc data
reordered PETSc initialization (now first part done in interface, setting of parameters done in numerics), removed unnecessary includes for PETSc still not running with gfortran, use at own risk!
2012-11-06 21:30:51 +05:30
DM, private :: da
numerics got some protected statements and is now reading in the new keywords for the solver selection in small letters. spectral solver got cut back facilities + improved output to screen
2012-10-02 20:56:56 +05:30
SNES, private :: snes
reordered PETSc initialization (now first part done in interface, setting of parameters done in numerics), removed unnecessary includes for PETSc still not running with gfortran, use at own risk!
2012-11-06 21:30:51 +05:30
Vec, private :: solution_vec
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
!--------------------------------------------------------------------------------------------------
! common pointwise data
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
real(pReal), private, dimension(:,:,:,:,:), allocatable :: &
F_lastInc, & !< field of previous compatible deformation gradients
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
F_lastInc2, & !< field of 2nd previous compatible deformation gradients
F_lambda_lastInc, & !< field of previous incompatible deformation gradient
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
Fdot, & !< field of assumed rate of compatible deformation gradient
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
F_lambdaDot !< field of assumed rate of incopatible deformation gradient
complex(pReal),private, dimension(:,:,:,:,:), allocatable :: inertiaField_fourier
fixed wrong temperature when using spectral solver
2013-01-09 23:38:08 +05:30
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
!--------------------------------------------------------------------------------------------------
! stress, stiffness and compliance average etc.
numerics got some protected statements and is now reading in the new keywords for the solver selection in small letters. spectral solver got cut back facilities + improved output to screen
2012-10-02 20:56:56 +05:30
real(pReal), private, dimension(3,3) :: &
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
F_aimDot, & !< assumed rate of average deformation gradient
F_aim = math_I3, & !< current prescribed deformation gradient
F_aim_lastInc = math_I3, & !< previous average deformation gradient
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
P_av = 0.0_pReal, & !< average 1st Piola--Kirchhoff stress
P_avLastEval = 0.0_pReal !< average 1st Piola--Kirchhoff stress last call of CPFEM_general
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
character(len=1024), private :: incInfo !< time and increment information
numerics got some protected statements and is now reading in the new keywords for the solver selection in small letters. spectral solver got cut back facilities + improved output to screen
2012-10-02 20:56:56 +05:30
real(pReal), private, dimension(3,3,3,3) :: &
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
C_volAvg = 0.0_pReal, & !< current volume average stiffness
C_volAvgLastInc = 0.0_pReal, & !< previous volume average stiffness
C_minMaxAvg = 0.0_pReal, & !< current (min+max)/2 stiffness
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
S = 0.0_pReal, & !< current compliance (filled up with zeros)
C_scale = 0.0_pReal, &
numerics got some protected statements and is now reading in the new keywords for the solver selection in small letters. spectral solver got cut back facilities + improved output to screen
2012-10-02 20:56:56 +05:30
S_scale = 0.0_pReal
next round in modularization
2012-07-25 19:31:39 +05:30
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
real(pReal), private :: &
err_stress, & !< deviation from stress BC
err_f, & !< difference between compatible and incompatible deformation gradient
err_p !< difference of stress resulting from compatible and incompatible F
logical, private :: ForwardData
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
integer(pInt), private :: &
totalIter = 0_pInt !< total iteration in current increment
added "public" (and where possible "protected") statements to functions/subroutines/variables/parameters where needed added -fmodule-private to compile options for gfortran, i.e. gfortran won't compile anymore if an entity is not explicitly declared as public but used in another module to ensure stronger encapsulation
2013-03-28 13:10:30 +05:30
public :: &
AL_init, &
AL_solution, &
AL_destroy
fixed bug in coordinate calculation when no average F was given (scaled wrong by ncp_elems**2) other changes: just polishing + some more comments
2013-03-27 17:58:55 +05:30
external :: &
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
VecDestroy, &
DMDestroy, &
DMDACreate3D, &
DMCreateGlobalVector, &
DMDASetLocalFunction, &
PETScFinalize, &
SNESDestroy, &
SNESGetNumberFunctionEvals, &
SNESGetIterationNumber, &
SNESSolve, &
SNESSetDM, &
SNESGetConvergedReason, &
SNESSetConvergenceTest, &
SNESSetFromOptions, &
SNESCreate, &
MPI_Abort
fixed bug in coordinate calculation when no average F was given (scaled wrong by ncp_elems**2) other changes: just polishing + some more comments
2013-03-27 17:58:55 +05:30
numerics got some protected statements and is now reading in the new keywords for the solver selection in small letters. spectral solver got cut back facilities + improved output to screen
2012-10-02 20:56:56 +05:30
contains
fixed bug in coordinate calculation when no average F was given (scaled wrong by ncp_elems**2) other changes: just polishing + some more comments
2013-03-27 17:58:55 +05:30
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief allocates all neccessary fields and fills them with data, potentially from restart info
made doxygen working for MSC.Marc again, small changes on the todo-statements for doxygen, they don't work on single lines of codes but only on module/variables
2013-06-11 22:05:04 +05:30
!> @todo use sourced allocation, e.g. allocate(Fdot,source = F_lastInc)
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
!--------------------------------------------------------------------------------------------------
fixed wrong temperature when using spectral solver
2013-01-09 23:38:08 +05:30
subroutine AL_init(temperature)
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
use, intrinsic :: iso_fortran_env ! to get compiler_version and compiler_options (at least for gfortran >4.6 at the moment)
use IO, only: &
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
IO_intOut, &
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
IO_read_JobBinaryFile, &
reporting current time for most module inits
2013-02-25 22:04:59 +05:30
IO_write_JobBinaryFile, &
IO_timeStamp
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
use debug, only : &
debug_level, &
debug_spectral, &
debug_spectralRestart
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
use FEsolving, only: &
restartInc
use DAMASK_interface, only: &
getSolverJobName
use DAMASK_spectral_Utilities, only: &
Utilities_init, &
Utilities_constitutiveResponse, &
Utilities_updateGamma, &
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
grid, &
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
grid1Red, &
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
geomSize, &
wgt
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
use mesh, only: &
added use statements
2013-02-05 20:33:36 +05:30
mesh_ipCoordinates, &
mesh_deformedCoordsFFT
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
use math, only: &
math_invSym3333
next round in modularization
2012-07-25 19:31:39 +05:30
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
implicit none
fixed wrong temperature when using spectral solver
2013-01-09 23:38:08 +05:30
real(pReal), intent(inout) :: &
temperature
some small changes in including PETSc *.h90 files to prevent warnings. It seems that a bug in old gfortran 4.5.0 causes the problems, at least with 4.7.3 it seems to work after fixing the bugs in IO.f90
2012-11-07 18:41:41 +05:30
#include <finclude/petscdmda.h90>
#include <finclude/petscsnes.h90>
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
real(pReal), dimension(:,:,:,:,:), allocatable :: P
made restart working for all solvers and added corresponding test corrected error handling for FFT based function in case of odd numbers etc.
2012-12-14 23:00:22 +05:30
real(pReal), dimension(3,3) :: &
temp33_Real = 0.0_pReal
real(pReal), dimension(3,3,3,3) :: &
temp3333_Real = 0.0_pReal, &
temp3333_Real2 = 0.0_pReal
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
PetscErrorCode :: ierr
PetscObject :: dummy
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
PetscScalar, pointer, dimension(:,:,:,:) :: xx_psc, F, F_lambda
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
call Utilities_init()
write(6,'(/,a)') ' <<<+- DAMASK_spectral_solverAL init -+>>>'
added svn properties
2013-08-01 14:43:46 +05:30
write(6,'(a)') ' $Id$'
reporting current time for most module inits
2013-02-25 22:04:59 +05:30
write(6,'(a16,a)') ' Current time : ',IO_timeStamp()
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)
2012-08-06 18:13:05 +05:30
#include "compilation_info.f90"
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
allocate (P (3,3,grid(1),grid(2),grid(3)),source = 0.0_pReal)
!--------------------------------------------------------------------------------------------------
! allocate global fields
made doxygen working for MSC.Marc again, small changes on the todo-statements for doxygen, they don't work on single lines of codes but only on module/variables
2013-06-11 22:05:04 +05:30
allocate (F_lastInc (3,3,grid(1),grid(2),grid(3)),source = 0.0_pReal)
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
allocate (F_lastInc2 (3,3,grid(1),grid(2),grid(3)),source = 0.0_pReal)
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
allocate (Fdot (3,3,grid(1),grid(2),grid(3)),source = 0.0_pReal)
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
allocate (F_lambda_lastInc(3,3,grid(1),grid(2),grid(3)),source = 0.0_pReal)
allocate (F_lambdaDot (3,3,grid(1),grid(2),grid(3)),source = 0.0_pReal)
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
allocate (inertiaField_fourier (grid1Red,grid(2),grid(3),3,3),source = cmplx(0.0_pReal,0.0_pReal,pReal))
changed order of most arrays to fortran-fast, whole new solver still work in progress
2012-08-09 18:34:56 +05:30
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
!--------------------------------------------------------------------------------------------------
! PETSc Init
loadcase rotation now working for AL solver
2012-12-17 15:48:39 +05:30
call SNESCreate(PETSC_COMM_WORLD,snes,ierr); CHKERRQ(ierr)
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
call DMDACreate3d(PETSC_COMM_WORLD, &
DMDA_BOUNDARY_NONE, DMDA_BOUNDARY_NONE, DMDA_BOUNDARY_NONE, &
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
DMDA_STENCIL_BOX,grid(1),grid(2),grid(3),PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE, &
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
18,1,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,da,ierr)
CHKERRQ(ierr)
loadcase rotation now working for AL solver
2012-12-17 15:48:39 +05:30
call DMCreateGlobalVector(da,solution_vec,ierr); CHKERRQ(ierr)
call DMDASetLocalFunction(da,AL_formResidual,ierr); CHKERRQ(ierr)
call SNESSetDM(snes,da,ierr); CHKERRQ(ierr)
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
call SNESSetConvergenceTest(snes,AL_converged,dummy,PETSC_NULL_FUNCTION,ierr)
CHKERRQ(ierr)
minor changes on buggy latest commit (wrong use statement)
2013-01-09 03:24:25 +05:30
call SNESSetFromOptions(snes,ierr); CHKERRQ(ierr)
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
!--------------------------------------------------------------------------------------------------
! init fields
loadcase rotation now working for AL solver
2012-12-17 15:48:39 +05:30
call DMDAVecGetArrayF90(da,solution_vec,xx_psc,ierr); CHKERRQ(ierr) ! places pointer xx_psc on PETSc data
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
F => xx_psc(0:8,:,:,:)
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
F_lambda => xx_psc(9:17,:,:,:)
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
if (restartInc == 1_pInt) then ! no deformation (no restart)
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
F_lastInc = spread(spread(spread(math_I3,3,grid(1)),4,grid(2)),5,grid(3)) ! initialize to identity
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
F_lastInc2 = F_lastInc
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
F_lambda_lastInc = F_lastInc
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
F = reshape(F_lastInc,[9,grid(1),grid(2),grid(3)])
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
F_lambda = F
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
elseif (restartInc > 1_pInt) then
if (iand(debug_level(debug_spectral),debug_spectralRestart)/= 0) &
write(6,'(/,a,'//IO_intOut(restartInc-1_pInt)//',a)') &
'reading values of increment', restartInc - 1_pInt, 'from file'
made restart working for all solvers and added corresponding test corrected error handling for FFT based function in case of odd numbers etc.
2012-12-14 23:00:22 +05:30
flush(6)
call IO_read_jobBinaryFile(777,'F',&
trim(getSolverJobName()),size(F))
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
read (777,rec=1) F
close (777)
made restart working for all solvers and added corresponding test corrected error handling for FFT based function in case of odd numbers etc.
2012-12-14 23:00:22 +05:30
call IO_read_jobBinaryFile(777,'F_lastInc',&
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
trim(getSolverJobName()),size(F_lastInc))
read (777,rec=1) F_lastInc
close (777)
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
call IO_read_jobBinaryFile(777,'F_lastInc2',&
trim(getSolverJobName()),size(F_lastInc2))
read (777,rec=1) F_lastInc2
close (777)
made restart working for all solvers and added corresponding test corrected error handling for FFT based function in case of odd numbers etc.
2012-12-14 23:00:22 +05:30
F_aim = reshape(sum(sum(sum(F,dim=4),dim=3),dim=2) * wgt, [3,3]) ! average of F
F_aim_lastInc = sum(sum(sum(F_lastInc,dim=5),dim=4),dim=3) * wgt ! average of F_lastInc
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
call IO_read_jobBinaryFile(777,'F_lambda',&
trim(getSolverJobName()),size(F_lambda))
read (777,rec=1) F_lambda
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
close (777)
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
call IO_read_jobBinaryFile(777,'F_lambda_lastInc',&
trim(getSolverJobName()),size(F_lambda_lastInc))
read (777,rec=1) F_lambda_lastInc
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
close (777)
introduced error code for run_test.py
2013-01-24 00:03:46 +05:30
call IO_read_jobBinaryFile(777,'F_aimDot',trim(getSolverJobName()),size(f_aimDot))
read (777,rec=1) f_aimDot
made restart working for all solvers and added corresponding test corrected error handling for FFT based function in case of odd numbers etc.
2012-12-14 23:00:22 +05:30
close (777)
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
call IO_read_jobBinaryFile(777,'C_volAvg',trim(getSolverJobName()),size(C_volAvg))
read (777,rec=1) C_volAvg
made restart working for all solvers and added corresponding test corrected error handling for FFT based function in case of odd numbers etc.
2012-12-14 23:00:22 +05:30
close (777)
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
call IO_read_jobBinaryFile(777,'C_volAvgLastInc',trim(getSolverJobName()),size(C_volAvgLastInc))
read (777,rec=1) C_volAvgLastInc
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
close (777)
made restart working for all solvers and added corresponding test corrected error handling for FFT based function in case of odd numbers etc.
2012-12-14 23:00:22 +05:30
call IO_read_jobBinaryFile(777,'C_ref',trim(getSolverJobName()),size(temp3333_Real))
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
read (777,rec=1) C_minMaxAvg
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
close (777)
endif
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
mesh_ipCoordinates = reshape(mesh_deformedCoordsFFT(geomSize,reshape(&
F,[3,3,grid(1),grid(2),grid(3)])),[3,1,product(grid)])
changed reference stiffness to min max avg from volume avg. using terminally ill information in convergence check
2013-03-06 20:01:13 +05:30
call Utilities_constitutiveResponse(F,F,temperature,0.0_pReal,P,temp3333_Real,temp3333_Real2,&
made restart working for all solvers and added corresponding test corrected error handling for FFT based function in case of odd numbers etc.
2012-12-14 23:00:22 +05:30
temp33_Real,.false.,math_I3)
loadcase rotation now working for AL solver
2012-12-17 15:48:39 +05:30
call DMDAVecRestoreArrayF90(da,solution_vec,xx_psc,ierr); CHKERRQ(ierr)
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
!--------------------------------------------------------------------------------------------------
! reference stiffness
made restart working for all solvers and added corresponding test corrected error handling for FFT based function in case of odd numbers etc.
2012-12-14 23:00:22 +05:30
if (restartInc == 1_pInt) then ! use initial stiffness as reference stiffness
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
C_minMaxAvg = temp3333_Real2
C_volAvg = temp3333_Real
made restart working for all solvers and added corresponding test corrected error handling for FFT based function in case of odd numbers etc.
2012-12-14 23:00:22 +05:30
endif
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
changed reference stiffness to min max avg from volume avg. using terminally ill information in convergence check
2013-03-06 20:01:13 +05:30
call Utilities_updateGamma(temp3333_Real2,.True.)
C_scale = temp3333_Real2
S_scale = math_invSym3333(temp3333_Real2)
next round in modularization
2012-07-25 19:31:39 +05:30
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
end subroutine AL_init
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief solution for the AL scheme with internal iterations
!--------------------------------------------------------------------------------------------------
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
type(tSolutionState) function &
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
AL_solution(incInfoIn,guess,timeinc,timeinc_old,loadCaseTime,P_BC,F_BC,temperature_bc,rotation_BC,density)
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
use numerics, only: &
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
update_gamma
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
use math, only: &
math_mul33x33 ,&
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
math_mul3333xx33, &
changed reference stiffness to min max avg from volume avg. using terminally ill information in convergence check
2013-03-06 20:01:13 +05:30
math_rotate_backward33, &
fixed wrong error calculation in Polarisation and AL solvers, and strictened tolerances for tests
2013-08-02 19:25:44 +05:30
math_invSym3333, &
math_transpose33
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
use mesh, only: &
added use statements
2013-02-05 20:33:36 +05:30
mesh_ipCoordinates, &
mesh_deformedCoordsFFT
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
use IO, only: &
IO_write_JobBinaryFile
use DAMASK_spectral_Utilities, only: &
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
grid, &
geomSize, &
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
tBoundaryCondition, &
Utilities_forwardField, &
Utilities_calculateRate, &
Utilities_maskedCompliance, &
Utilities_updateGamma, &
cutBack
use FEsolving, only: &
restartWrite, &
terminallyIll
implicit none
some small changes in including PETSc *.h90 files to prevent warnings. It seems that a bug in old gfortran 4.5.0 causes the problems, at least with 4.7.3 it seems to work after fixing the bugs in IO.f90
2012-11-07 18:41:41 +05:30
#include <finclude/petscdmda.h90>
#include <finclude/petscsnes.h90>
introduced possibility so specify deformation gradient aim at end of load case, rate will be calculated using difference between start of load case and aim. Needed for cyclic loading. Use keyword "f" for this behavior, don't use it as short name for freq any more!
2013-05-13 15:14:23 +05:30
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
!--------------------------------------------------------------------------------------------------
! input data for solution
introduced possibility so specify deformation gradient aim at end of load case, rate will be calculated using difference between start of load case and aim. Needed for cyclic loading. Use keyword "f" for this behavior, don't use it as short name for freq any more!
2013-05-13 15:14:23 +05:30
real(pReal), intent(in) :: &
timeinc, & !< increment in time for current solution
timeinc_old, & !< increment in time of last increment
loadCaseTime, & !< remaining time of current load case
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
temperature_bc, &
density
loadcase rotation now working for AL solver
2012-12-17 15:48:39 +05:30
logical, intent(in) :: &
guess
type(tBoundaryCondition), intent(in) :: &
P_BC, &
F_BC
character(len=*), intent(in) :: &
incInfoIn
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
real(pReal), dimension(3,3), intent(in) :: rotation_BC
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
!--------------------------------------------------------------------------------------------------
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
! PETSc Data
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
PetscScalar, dimension(:,:,:,:), pointer :: xx_psc, F, F_lambda
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
PetscErrorCode :: ierr
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
SNESConvergedReason ::reason
more work on array ordering in new solver. (now working but not tested)
2012-08-10 22:31:58 +05:30
fixed wrong error calculation in Polarisation and AL solvers, and strictened tolerances for tests
2013-08-02 19:25:44 +05:30
integer(pInt) :: i, j, k
real(pReal), dimension(3,3) :: F_lambda33
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
incInfo = incInfoIn
made restart working for all solvers and added corresponding test corrected error handling for FFT based function in case of odd numbers etc.
2012-12-14 23:00:22 +05:30
call DMDAVecGetArrayF90(da,solution_vec,xx_psc,ierr)
F => xx_psc(0:8,:,:,:)
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
F_lambda => xx_psc(9:17,:,:,:)
made restart working for all solvers and added corresponding test corrected error handling for FFT based function in case of odd numbers etc.
2012-12-14 23:00:22 +05:30
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
!--------------------------------------------------------------------------------------------------
! restart information for spectral solver
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
if (restartWrite) then
corrected reporting of basic PETSc variant and unified reporting of all spectral solvers improved warning and error in IO, now able to report correctly elements up to 9 digits
2013-01-10 03:49:32 +05:30
write(6,'(/,a)') ' writing converged results for restart'
flush(6)
made restart working for all solvers and added corresponding test corrected error handling for FFT based function in case of odd numbers etc.
2012-12-14 23:00:22 +05:30
call IO_write_jobBinaryFile(777,'F',size(F)) ! writing deformation gradient field to file
write (777,rec=1) F
close (777)
call IO_write_jobBinaryFile(777,'F_lastInc',size(F_lastInc)) ! writing F_lastInc field to file
write (777,rec=1) F_lastInc
close (777)
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
call IO_write_jobBinaryFile(777,'F_lambda',size(F_lambda)) ! writing deformation gradient field to file
write (777,rec=1) F_lambda
made restart working for all solvers and added corresponding test corrected error handling for FFT based function in case of odd numbers etc.
2012-12-14 23:00:22 +05:30
close (777)
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
call IO_write_jobBinaryFile(777,'F_lambda_lastInc',size(F_lambda_lastInc)) ! writing F_lastInc field to file
write (777,rec=1) F_lambda_lastInc
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
close (777)
added doxygen documentation to material.f90 and marked read-only quantities as protected where possible, removed substituted "call flush" by "flush"
2013-01-18 17:00:52 +05:30
call IO_write_jobBinaryFile(777,'F_aimDot',size(F_aimDot))
write (777,rec=1) F_aimDot
close(777)
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
call IO_write_jobBinaryFile(777,'C_volAvg',size(C_volAvg))
write (777,rec=1) C_volAvg
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
close(777)
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
call IO_write_jobBinaryFile(777,'C_volAvgLastInc',size(C_volAvgLastInc))
write (777,rec=1) C_volAvgLastInc
made restart working for all solvers and added corresponding test corrected error handling for FFT based function in case of odd numbers etc.
2012-12-14 23:00:22 +05:30
close(777)
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
endif
AL_solution%converged =.false.
introduced possibility so specify deformation gradient aim at end of load case, rate will be calculated using difference between start of load case and aim. Needed for cyclic loading. Use keyword "f" for this behavior, don't use it as short name for freq any more!
2013-05-13 15:14:23 +05:30
if (cutBack) then
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
F_aim = F_aim_lastInc
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
F_lambda= reshape(F_lambda_lastInc,[9,grid(1),grid(2),grid(3)])
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
F = reshape(F_lastInc, [9,grid(1),grid(2),grid(3)])
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
C_volAvg = C_volAvgLastInc
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
else
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
C_volAvgLastInc = C_volAvg
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
!--------------------------------------------------------------------------------------------------
numerics got some protected statements and is now reading in the new keywords for the solver selection in small letters. spectral solver got cut back facilities + improved output to screen
2012-10-02 20:56:56 +05:30
! calculate rate for aim
if (F_BC%myType=='l') then ! calculate f_aimDot from given L and current F
f_aimDot = F_BC%maskFloat * math_mul33x33(F_BC%values, F_aim)
introduced possibility so specify deformation gradient aim at end of load case, rate will be calculated using difference between start of load case and aim. Needed for cyclic loading. Use keyword "f" for this behavior, don't use it as short name for freq any more!
2013-05-13 15:14:23 +05:30
elseif(F_BC%myType=='fdot') then ! f_aimDot is prescribed
numerics got some protected statements and is now reading in the new keywords for the solver selection in small letters. spectral solver got cut back facilities + improved output to screen
2012-10-02 20:56:56 +05:30
f_aimDot = F_BC%maskFloat * F_BC%values
introduced possibility so specify deformation gradient aim at end of load case, rate will be calculated using difference between start of load case and aim. Needed for cyclic loading. Use keyword "f" for this behavior, don't use it as short name for freq any more!
2013-05-13 15:14:23 +05:30
elseif(F_BC%myType=='f') then ! aim at end of load case is prescribed
f_aimDot = F_BC%maskFloat * (F_BC%values -F_aim)/loadCaseTime
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
endif
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
if (guess) f_aimDot = f_aimDot + P_BC%maskFloat * (F_aim - F_aim_lastInc)/timeinc_old
numerics got some protected statements and is now reading in the new keywords for the solver selection in small letters. spectral solver got cut back facilities + improved output to screen
2012-10-02 20:56:56 +05:30
F_aim_lastInc = F_aim
!--------------------------------------------------------------------------------------------------
! update coordinates and rate and forward last inc
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
mesh_ipCoordinates = reshape(mesh_deformedCoordsFFT(geomSize,reshape(&
F,[3,3,grid(1),grid(2),grid(3)])),[3,1,product(grid)])
numerics got some protected statements and is now reading in the new keywords for the solver selection in small letters. spectral solver got cut back facilities + improved output to screen
2012-10-02 20:56:56 +05:30
Fdot = Utilities_calculateRate(math_rotate_backward33(f_aimDot,rotation_BC), &
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
timeinc_old,guess,F_lastInc,reshape(F,[3,3,grid(1),grid(2),grid(3)]))
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
F_lambdaDot = Utilities_calculateRate(math_rotate_backward33(f_aimDot,rotation_BC), &
timeinc_old,guess,F_lambda_lastInc,reshape(F_lambda,[3,3,grid(1),grid(2),grid(3)]))
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
F_lastInc2 = F_lastInc
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
F_lastInc = reshape(F, [3,3,grid(1),grid(2),grid(3)])
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
F_lambda_lastInc = reshape(F_lambda,[3,3,grid(1),grid(2),grid(3)])
numerics got some protected statements and is now reading in the new keywords for the solver selection in small letters. spectral solver got cut back facilities + improved output to screen
2012-10-02 20:56:56 +05:30
endif
F_aim = F_aim + f_aimDot * timeinc
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
!--------------------------------------------------------------------------------------------------
loadcase rotation now working for AL solver
2012-12-17 15:48:39 +05:30
! update local deformation gradient
added debugging possibility for MSC.Marc, rename parameters to CAPITALS
2013-03-31 18:36:49 +05:30
F = reshape(Utilities_forwardField(timeinc,F_lastInc,Fdot, & ! ensure that it matches rotated F_aim
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
math_rotate_backward33(F_aim,rotation_BC)),[9,grid(1),grid(2),grid(3)])
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
F_lambda = reshape(Utilities_forwardField(timeinc,F_lambda_lastInc,F_lambdadot), [9,grid(1),grid(2),grid(3)]) ! does not have any average value as boundary condition
fixed wrong error calculation in Polarisation and AL solvers, and strictened tolerances for tests
2013-08-02 19:25:44 +05:30
if (.not. guess) then ! large strain forwarding
do k = 1_pInt, grid(3); do j = 1_pInt, grid(2); do i = 1_pInt, grid(1)
F_lambda33 = reshape(F_lambda(:,i,j,k),[3,3])
F_lambda33 = math_mul3333xx33(S_scale,math_mul33x33(F_lambda33, &
math_mul3333xx33(C_scale,&
math_mul33x33(math_transpose33(F_lambda33),&
F_lambda33) -math_I3))*0.5_pReal)&
+ math_I3
F_lambda(:,i,j,k) = reshape(F_lambda33,[9])
enddo; enddo; enddo
endif
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
call DMDAVecRestoreArrayF90(da,solution_vec,xx_psc,ierr)
CHKERRQ(ierr)
numerics got some protected statements and is now reading in the new keywords for the solver selection in small letters. spectral solver got cut back facilities + improved output to screen
2012-10-02 20:56:56 +05:30
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
!--------------------------------------------------------------------------------------------------
! update stiffness (and gamma operator)
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
S = Utilities_maskedCompliance(rotation_BC,P_BC%maskLogical,C_volAvg)
changed reference stiffness to min max avg from volume avg. using terminally ill information in convergence check
2013-03-06 20:01:13 +05:30
if (update_gamma) then
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
call Utilities_updateGamma(C_minMaxAvg,restartWrite)
C_scale = C_minMaxAvg
S_scale = math_invSym3333(C_minMaxAvg)
changed reference stiffness to min max avg from volume avg. using terminally ill information in convergence check
2013-03-06 20:01:13 +05:30
endif
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
ForwardData = .True.
mask_stress = P_BC%maskFloat
params%P_BC = P_BC%values
params%rotation_BC = rotation_BC
params%timeinc = timeinc
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
params%timeincOld = timeinc_old
fixed wrong temperature when using spectral solver
2013-01-09 23:38:08 +05:30
params%temperature = temperature_bc
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
params%density = density
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
!--------------------------------------------------------------------------------------------------
! solve BVP
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
call SNESSolve(snes,PETSC_NULL_OBJECT,solution_vec,ierr)
CHKERRQ(ierr)
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
!--------------------------------------------------------------------------------------------------
! check convergence
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
call SNESGetConvergedReason(snes,reason,ierr)
CHKERRQ(ierr)
AL_solution%termIll = terminallyIll
terminallyIll = .false.
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
AL_solution%converged = .true.
if (reason < 1 ) AL_solution%converged = .false.
AL_solution%iterationsNeeded = totalIter
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
end function AL_solution
numerics got some protected statements and is now reading in the new keywords for the solver selection in small letters. spectral solver got cut back facilities + improved output to screen
2012-10-02 20:56:56 +05:30
new version of modular solver structure
2012-08-03 14:55:48 +05:30
moved petsc option string to numerics and added new numerics parameters to config file
2012-08-06 22:57:53 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief forms the AL residual vector
!--------------------------------------------------------------------------------------------------
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
subroutine AL_formResidual(in,x_scal,f_scal,dummy,ierr)
use numerics, only: &
itmax, &
restructured AL solver to resemble polarization scheme and introduced alpha and beta parameters in numerics to determine what kind of solution you want. alpha = beta = 1 results in the AL scheme. alpha = beta = 2 results in the accelerated scheme. safe values are anything in between. changed the solver to use stress calculations on the campatible field which seems to result in faster convergence.
2013-02-28 23:05:02 +05:30
itmin, &
polarAlpha, &
fixed wrong error calculation in Polarisation and AL solvers, and strictened tolerances for tests
2013-08-02 19:25:44 +05:30
polarBeta
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
use IO, only: &
IO_intOut
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
use math, only: &
math_rotate_backward33, &
math_transpose33, &
renamed F_lambda to F_tau for clarity. using stress tangent to precondition equations
2013-03-04 15:19:40 +05:30
math_mul3333xx33, &
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
math_invSym3333, &
PI
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
use DAMASK_spectral_Utilities, only: &
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
grid, &
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
geomSize, &
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
wgt, &
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
field_real, &
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
field_fourier, &
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
Utilities_FFTforward, &
Utilities_fourierConvolution, &
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
Utilities_inverseLaplace, &
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
Utilities_FFTbackward, &
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
Utilities_constitutiveResponse
use debug, only: &
debug_level, &
debug_spectral, &
debug_spectralRotation
renamed F_lambda to F_tau for clarity. using stress tangent to precondition equations
2013-03-04 15:19:40 +05:30
use homogenization, only: &
materialpoint_dPdF
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
implicit none
corrected PETSc macro expansion
2013-01-10 21:06:55 +05:30
!--------------------------------------------------------------------------------------------------
! strange syntax in the next line because otherwise macros expand beyond 132 character limit
DMDALocalInfo, dimension(&
DMDA_LOCAL_INFO_SIZE) :: &
added some warning and explicit size of arrays
2013-01-10 19:03:43 +05:30
in
corrected PETSc macro expansion
2013-01-10 21:06:55 +05:30
PetscScalar, target, dimension(3,3,2, &
XG_RANGE,YG_RANGE,ZG_RANGE) :: &
added some warning and explicit size of arrays
2013-01-10 19:03:43 +05:30
x_scal
corrected PETSc macro expansion
2013-01-10 21:06:55 +05:30
PetscScalar, target, dimension(3,3,2, &
X_RANGE,Y_RANGE,Z_RANGE) :: &
added some warning and explicit size of arrays
2013-01-10 19:03:43 +05:30
f_scal
corrected PETSc macro expansion
2013-01-10 21:06:55 +05:30
PetscScalar, pointer, dimension(:,:,:,:,:) :: &
added some warning and explicit size of arrays
2013-01-10 19:03:43 +05:30
F, &
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
F_lambda, &
corrected PETSc macro expansion
2013-01-10 21:06:55 +05:30
residual_F, &
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
residual_F_lambda
added some warning and explicit size of arrays
2013-01-10 19:03:43 +05:30
PetscInt :: &
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
PETScIter, &
added some warning and explicit size of arrays
2013-01-10 19:03:43 +05:30
nfuncs
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
PetscObject :: dummy
PetscErrorCode :: ierr
forgot some integers
2013-01-11 00:20:14 +05:30
integer(pInt) :: &
fixed wrong error calculation in Polarisation and AL solvers, and strictened tolerances for tests
2013-08-02 19:25:44 +05:30
i, j, k, e
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
real(pReal) :: correctionFactor
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
added debugging possibility for MSC.Marc, rename parameters to CAPITALS
2013-03-31 18:36:49 +05:30
F => x_scal(1:3,1:3,1,&
corrected PETSc macro expansion
2013-01-10 21:06:55 +05:30
XG_RANGE,YG_RANGE,ZG_RANGE)
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
F_lambda => x_scal(1:3,1:3,2,&
corrected PETSc macro expansion
2013-01-10 21:06:55 +05:30
XG_RANGE,YG_RANGE,ZG_RANGE)
added debugging possibility for MSC.Marc, rename parameters to CAPITALS
2013-03-31 18:36:49 +05:30
residual_F => f_scal(1:3,1:3,1,&
corrected PETSc macro expansion
2013-01-10 21:06:55 +05:30
X_RANGE,Y_RANGE,Z_RANGE)
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
residual_F_lambda => f_scal(1:3,1:3,2,&
corrected PETSc macro expansion
2013-01-10 21:06:55 +05:30
X_RANGE,Y_RANGE,Z_RANGE)
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
corrected reporting of basic PETSc variant and unified reporting of all spectral solvers improved warning and error in IO, now able to report correctly elements up to 9 digits
2013-01-10 03:49:32 +05:30
call SNESGetNumberFunctionEvals(snes,nfuncs,ierr); CHKERRQ(ierr)
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
call SNESGetIterationNumber(snes,PETScIter,ierr); CHKERRQ(ierr)
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
if(nfuncs== 0 .and. PETScIter == 0) totalIter = -1_pInt ! new increment
if (totalIter <= PETScIter) then ! new iteration
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
!--------------------------------------------------------------------------------------------------
! report begin of new iteration
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
totalIter = totalIter + 1_pInt
corrected reporting of basic PETSc variant and unified reporting of all spectral solvers improved warning and error in IO, now able to report correctly elements up to 9 digits
2013-01-10 03:49:32 +05:30
write(6,'(1x,a,3(a,'//IO_intOut(itmax)//'))') trim(incInfo), &
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
' @ Iteration ', itmin, '≤',totalIter, '≤', itmax
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
if (iand(debug_level(debug_spectral),debug_spectralRotation) /= 0) &
corrected reporting of basic PETSc variant and unified reporting of all spectral solvers improved warning and error in IO, now able to report correctly elements up to 9 digits
2013-01-10 03:49:32 +05:30
write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') ' deformation gradient aim (lab) =', &
loadcase rotation now working for AL solver
2012-12-17 15:48:39 +05:30
math_transpose33(math_rotate_backward33(F_aim,params%rotation_BC))
corrected reporting of basic PETSc variant and unified reporting of all spectral solvers improved warning and error in IO, now able to report correctly elements up to 9 digits
2013-01-10 03:49:32 +05:30
write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') ' deformation gradient aim =', &
loadcase rotation now working for AL solver
2012-12-17 15:48:39 +05:30
math_transpose33(F_aim)
corrected reporting of basic PETSc variant and unified reporting of all spectral solvers improved warning and error in IO, now able to report correctly elements up to 9 digits
2013-01-10 03:49:32 +05:30
flush(6)
numerics got some protected statements and is now reading in the new keywords for the solver selection in small letters. spectral solver got cut back facilities + improved output to screen
2012-10-02 20:56:56 +05:30
endif
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
!--------------------------------------------------------------------------------------------------
! evaluate inertia
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
dynamic: if (params%density > 0.0_pReal) then
residual_F = ((F - F_lastInc)/params%timeinc - (F_lastInc - F_lastInc2)/params%timeincOld)/&
((params%timeinc + params%timeincOld)/2.0_pReal)
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
residual_F = params%density*product(geomSize/grid)*residual_F
field_real = 0.0_pReal
field_real(1:grid(1),1:grid(2),1:grid(3),1:3,1:3) = reshape(residual_F,[grid(1),grid(2),grid(3),3,3],&
order=[4,5,1,2,3]) ! field real has a different order
call Utilities_FFTforward()
call Utilities_inverseLaplace()
inertiaField_fourier = field_fourier
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
else dynamic
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
inertiaField_fourier = cmplx(0.0_pReal,0.0_pReal,pReal)
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
endif dynamic
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
!--------------------------------------------------------------------------------------------------
!
field_real = 0.0_pReal
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
do k = 1_pInt, grid(3); do j = 1_pInt, grid(2); do i = 1_pInt, grid(1)
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
field_real(i,j,k,1:3,1:3) = math_mul3333xx33(C_scale,polarBeta*F(1:3,1:3,i,j,k) - &
polarAlpha*F_lambda(1:3,1:3,i,j,k))
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
enddo; enddo; enddo
!--------------------------------------------------------------------------------------------------
! doing convolution in Fourier space
call Utilities_FFTforward()
introduced dynamics for high strain rate simulations. so far only works for materials with homogeneous density. the density is read in as an optional argument in the load file using the keyword 'den' or 'density'. if density > 0 then inertial terms are added to the stress equilibrium. an implicit euler time discretisation is used to calculate the inertial terms. dynamic problems are length scale dependent so make sure the geometry size is meaningful (1 m^3 or N m^3 for example is too big for a RVE!). not specifying a density (density = 0) will perform a quasi static simulation as before.
2013-07-26 21:55:37 +05:30
field_fourier = field_fourier + polarAlpha*inertiaField_fourier
restructured AL solver to resemble polarization scheme and introduced alpha and beta parameters in numerics to determine what kind of solution you want. alpha = beta = 1 results in the AL scheme. alpha = beta = 2 results in the accelerated scheme. safe values are anything in between. changed the solver to use stress calculations on the campatible field which seems to result in faster convergence.
2013-02-28 23:05:02 +05:30
call Utilities_fourierConvolution(math_rotate_backward33(polarBeta*F_aim,params%rotation_BC))
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
call Utilities_FFTbackward()
!--------------------------------------------------------------------------------------------------
! constructing residual
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
residual_F_lambda = polarBeta*F - reshape(field_real(1:grid(1),1:grid(2),1:grid(3),1:3,1:3),&
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
[3,3,grid(1),grid(2),grid(3)],order=[3,4,5,1,2])
restructured AL solver to resemble polarization scheme and introduced alpha and beta parameters in numerics to determine what kind of solution you want. alpha = beta = 1 results in the AL scheme. alpha = beta = 2 results in the accelerated scheme. safe values are anything in between. changed the solver to use stress calculations on the campatible field which seems to result in faster convergence.
2013-02-28 23:05:02 +05:30
!--------------------------------------------------------------------------------------------------
! evaluate constitutive response
fixed bug in AL, weakened correction slightly
2013-07-31 02:34:41 +05:30
P_avLastEval = P_av
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
call Utilities_constitutiveResponse(F_lastInc,F - residual_F_lambda/polarBeta,params%temperature,params%timeinc, &
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
residual_F,C_volAvg,C_minMaxAvg,P_av,ForwardData,params%rotation_BC)
restructured AL solver to resemble polarization scheme and introduced alpha and beta parameters in numerics to determine what kind of solution you want. alpha = beta = 1 results in the AL scheme. alpha = beta = 2 results in the accelerated scheme. safe values are anything in between. changed the solver to use stress calculations on the campatible field which seems to result in faster convergence.
2013-02-28 23:05:02 +05:30
ForwardData = .False.
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
!--------------------------------------------------------------------------------------------------
! stress BC handling
write(6,'(/,a)') ' ... correcting F to fullfill stress BC ....................................'
fixed bug in AL, weakened correction slightly
2013-07-31 02:34:41 +05:30
correctionFactor = (cos((1.0-10000.0_pReal**(-sum((P_av-P_avLastEval)**2.0_pReal)/& ! only correct when averages stress of last two calls doesn't strongly deviate
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
sum(P_av**2.0_pReal)))*PI)+1.0)/2.0_pReal
fixed bug in AL, weakened correction slightly
2013-07-31 02:34:41 +05:30
write(6,'(/,a,f10.4)') ' stress BC correction factor = ', correctionFactor
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
F_aim = F_aim - correctionFactor *math_mul3333xx33(S, ((P_av - params%P_BC))) ! S = 0.0 for no bc
err_stress = maxval(abs(mask_stress * (P_av - params%P_BC))) ! mask = 0.0 for no bc
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
restructured AL solver to resemble polarization scheme and introduced alpha and beta parameters in numerics to determine what kind of solution you want. alpha = beta = 1 results in the AL scheme. alpha = beta = 2 results in the accelerated scheme. safe values are anything in between. changed the solver to use stress calculations on the campatible field which seems to result in faster convergence.
2013-02-28 23:05:02 +05:30
!--------------------------------------------------------------------------------------------------
! constructing residual
fixed wrong error calculation in Polarisation and AL solvers, and strictened tolerances for tests
2013-08-02 19:25:44 +05:30
e = 0_pInt
renamed F_lambda to F_tau for clarity. using stress tangent to precondition equations
2013-03-04 15:19:40 +05:30
err_p = 0.0_pReal
moved public data res,size and homog from mesh to DAMASK_spectral_utilities (as grid and geomSize)
2013-05-08 21:22:29 +05:30
do k = 1_pInt, grid(3); do j = 1_pInt, grid(2); do i = 1_pInt, grid(1)
fixed wrong error calculation in Polarisation and AL solvers, and strictened tolerances for tests
2013-08-02 19:25:44 +05:30
e = e + 1_pInt
residual_F(1:3,1:3,i,j,k) = math_mul3333xx33(math_invSym3333(materialpoint_dPdF(:,:,:,:,1,e) + C_scale), &
residual_F(1:3,1:3,i,j,k) - &
math_mul3333xx33(C_scale,F_lambda(1:3,1:3,i,j,k) - math_I3)) &
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
+ residual_F_lambda(1:3,1:3,i,j,k)
corrected invalid error for mismatch_p
2013-08-02 17:18:24 +05:30
err_p = err_p + sum((math_mul3333xx33(C_scale,residual_F(1:3,1:3,i,j,k) - residual_F_lambda(1:3,1:3,i,j,k)))**2.0_pReal)
restructured AL solver to resemble polarization scheme and introduced alpha and beta parameters in numerics to determine what kind of solution you want. alpha = beta = 1 results in the AL scheme. alpha = beta = 2 results in the accelerated scheme. safe values are anything in between. changed the solver to use stress calculations on the campatible field which seems to result in faster convergence.
2013-02-28 23:05:02 +05:30
enddo; enddo; enddo
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
!--------------------------------------------------------------------------------------------------
! calculating errors
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
err_f = wgt*sqrt(sum(residual_F_lambda**2.0_pReal))/polarBeta
renamed F_lambda to F_tau for clarity. using stress tangent to precondition equations
2013-03-04 15:19:40 +05:30
err_p = wgt*sqrt(err_p)
changed order of most arrays to fortran-fast, whole new solver still work in progress
2012-08-09 18:34:56 +05:30
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
end subroutine AL_formResidual
updates for AL version of spectral solver
2012-08-06 14:23:12 +05:30
moved petsc option string to numerics and added new numerics parameters to config file
2012-08-06 22:57:53 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief convergence check
!--------------------------------------------------------------------------------------------------
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
subroutine AL_converged(snes_local,PETScIter,xnorm,snorm,fnorm,reason,dummy,ierr)
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
use numerics, only: &
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
itmax, &
itmin, &
set tolerances to useful values and renamed them. allowed "keyword = value" in numerics.config (simply removing "=" from line)
2013-07-31 20:33:38 +05:30
err_f_tolabs, &
err_p_tolabs, &
err_f_p_tolrel, &
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
err_stress_tolabs, &
err_stress_tolrel
changed reference stiffness to min max avg from volume avg. using terminally ill information in convergence check
2013-03-06 20:01:13 +05:30
use FEsolving, only: &
terminallyIll
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
corrected reporting of basic PETSc variant and unified reporting of all spectral solvers improved warning and error in IO, now able to report correctly elements up to 9 digits
2013-01-10 03:49:32 +05:30
implicit none
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
SNES :: snes_local
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
PetscInt :: PETScIter
corrected reporting of basic PETSc variant and unified reporting of all spectral solvers improved warning and error in IO, now able to report correctly elements up to 9 digits
2013-01-10 03:49:32 +05:30
PetscReal :: &
xnorm, &
snorm, &
fnorm
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
SNESConvergedReason :: reason
PetscObject :: dummy
PetscErrorCode ::ierr
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
real(pReal) :: &
mismatch_f_tol, &
mismatch_p_tol, &
stressBC_tol
set tolerances to useful values and renamed them. allowed "keyword = value" in numerics.config (simply removing "=" from line)
2013-07-31 20:33:38 +05:30
mismatch_f_tol = max(maxval(abs(F_aim-math_I3))*err_f_p_tolrel,err_f_tolabs)
mismatch_p_tol = max(maxval(abs(P_av)) *err_f_p_tolrel,err_p_tolabs)
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
stressBC_tol = max(maxval(abs(P_av)) *err_stress_tolrel,err_stress_tolabs)
decoupled compatibility/equilibrium calculation from stress BC correction. needs more iteration for mildly contrasted materials, but there the basic scheme is better suited anyway. but now converges better for highly contrasted VEs
2013-07-23 23:12:15 +05:30
write(6,'(1/,a)') ' ... reporting .............................................................'
corrected reporting of basic PETSc variant and unified reporting of all spectral solvers improved warning and error in IO, now able to report correctly elements up to 9 digits
2013-01-10 03:49:32 +05:30
write(6,'(/,a,f8.2,a,es11.5,a,es11.4,a)') ' mismatch F = ', &
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
err_f/mismatch_f_tol, &
' (',err_f,' -, tol =',mismatch_f_tol,')'
corrected reporting of basic PETSc variant and unified reporting of all spectral solvers improved warning and error in IO, now able to report correctly elements up to 9 digits
2013-01-10 03:49:32 +05:30
write(6,'(a,f8.2,a,es11.5,a,es11.4,a)') ' mismatch P = ', &
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
err_p/mismatch_p_tol, &
' (',err_p,' -, tol =',mismatch_p_tol,')'
write(6,'(a,f8.2,a,es11.5,a,es11.4,a)') ' error stress BC = ', &
err_stress/stressBC_tol, ' (',err_stress, ' Pa, tol =',stressBC_tol,')'
write(6,'(/,a)') ' ==========================================================================='
flush(6)
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
converged: if ((totalIter >= itmin .and. &
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
all([ err_f/mismatch_f_tol, &
err_p/mismatch_p_tol, &
err_stress/stressBC_tol] < 1.0_pReal)) &
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
.or. terminallyIll) then
reason = 1
changes on convergence tolerances of AL and Polarisation, switched back to immediate correction of stress bc but only when last two average stresses are close to each other (cosine decay)
2013-07-30 21:02:55 +05:30
elseif (totalIter >= itmax) then converged
added Polarisation scheme as solver
2013-07-24 18:36:16 +05:30
reason = -1
else converged
reason = 0
endif converged
corrected reporting of basic PETSc variant and unified reporting of all spectral solvers improved warning and error in IO, now able to report correctly elements up to 9 digits
2013-01-10 03:49:32 +05:30
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
end subroutine AL_converged
new version of modular solver structure
2012-08-03 14:55:48 +05:30
moved petsc option string to numerics and added new numerics parameters to config file
2012-08-06 22:57:53 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief destroy routine
!--------------------------------------------------------------------------------------------------
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
subroutine AL_destroy()
use DAMASK_spectral_Utilities, only: &
Utilities_destroy
corrected reporting of basic PETSc variant and unified reporting of all spectral solvers improved warning and error in IO, now able to report correctly elements up to 9 digits
2013-01-10 03:49:32 +05:30
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
implicit none
PetscErrorCode :: ierr
corrected reporting of basic PETSc variant and unified reporting of all spectral solvers improved warning and error in IO, now able to report correctly elements up to 9 digits
2013-01-10 03:49:32 +05:30
call VecDestroy(solution_vec,ierr); CHKERRQ(ierr)
call SNESDestroy(snes,ierr); CHKERRQ(ierr)
call DMDestroy(da,ierr); CHKERRQ(ierr)
call PetscFinalize(ierr); CHKERRQ(ierr)
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
call Utilities_destroy()
more modularization….
2012-07-26 19:28:47 +05:30
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
end subroutine AL_destroy
more modularization….
2012-07-26 19:28:47 +05:30
added comments and structured the code. Temperature is not longer stored for each point but only to simplify restart behavior as long as it is not fully supported by the constitutive models
2012-11-12 19:44:39 +05:30
end module DAMASK_spectral_SolverAL