Programming
FEAP is designed to allow users to modify the code for their owm purposes. If at all possible, users should always use the user hooks provided in the user directory. These allow users to make custom elements, material models, new FEAP commands, et cetera. The programmers manual provides an introduction to programming in FEAP.
To use the user hooks, it is recommended to make a copy of the relevant stub file in a separate project directory. Then add the full path name to this edited copy to main/makefile on the OBJECTS line so that it will be built into FEAP the next time you compile it. If you are using windows, add the copy to your executable project. Even though there is a version of the stub file in your FEAP archive file, the compiler should use your edited copy instead of the one in the archive.
It is also recommended that you obtain a good Fortran reference book, e.g. Modern Fortran Explained by Metcalf, Reid, and Cohen.
Programmers Manual
The programmers manual can be found on the FEAP project site http://projects.ce.berkeley.edu/feap.
User Elements
User Interface Elements
User Memory Allocation
User Block Generation
User Body Forces
User Materials
The method to implement a user material model in FEAP has been described in Section 4.2 of FEAP Programmer Manual. Here we would like to clarify a few things in the user material model with history variables.
The time-dependent history variables
The time-dependent variables are meant to be used for quantities that change in "TIME". That is if you have a solution that advances from time to the history variables for are stored in hr(nh1)
and those for are stored in hr(nh2)
. These can be state variables in incremental stress-strain relations. In your solution algorithm, you must use the solution commands TIME
and DT
. The structure of your user material model should look like this:
1. In your umati*.f routine, set n1
equals the number of history variabales per Gauss point.
2. Read the history variables from the previous globally converged time step from hn(*)
3. Compute new values of history variables
4. Write your newly computed history values to h1(*)
for next time step . Do not write to hn(*)
!
After that, FEAP will carry out a global Newton iteration. If the global Newton iteration is not converged, then the values you wrote to h1(*)
will be discarded and your routine will be called again. You will repeat again from step 2 until it is converged. Then, the values in h1(*)
will be copied over into hn(*)
for use in the next time step when you execute the next TIME
command.
The time-independent history variables
If the history variable does not depend in time then we consider them to be independent of "TIME", so you could use the time-independent history variables as it only stores one copy of each value instead of two in the previous case. This type of history variable is more like a status variable depending on deformation. Here is the procedure to implement such kind of history variables,
1. In your umati*.f routine, set n3
equals the number of history variables per element.
2. Read the history variables from the previous globally converged time step from hr(nh3+i)
where i
runs from 0 to the number of history variables minus one and nh3
points to the location of first history variable for the current Gauss point.
3. Compute new values of history variables.
4. Write your newly computed history values to hr(*)
if needed.
Initialization of history variables in element or materials routines
The history variables are initialized to zero by default. But if you want to initialize any of the history variables to non-zero values you can do this when ISW .eq. 14
. For example, if we have m
history variables per gauss point, then you can initialize the time-dependent history variables to 1.0d0
in your umodelf.f or umodelf.f subroutine as follows,
c... INITIALIZE HISTORY VARIABLES IF (ISW .EQ. 14) THEN do i = 1, m hn(i) = 1.0d0 h1(i) = 1.0d0 end do END IF
For time-independent history variables
c... INITIALIZE HISTORY VARIABLES IF (ISW .EQ. 14) THEN do i = 0, (m - 1) hr(nh3 + i ) = 1.0d0 end do END IF
This segment of code will be executed for each Gauss point of each element.
(Note that this page is still under development)