Here is how I would do it:

(1) For each integration point, when the VUMAT is called for the first time, read the euler angles (they are stored in the PROPS vector) and calculate the initial rotation matrix (at the beginning of the time increment) that would allow to rotate the single crystal properties in the appropriate coordinate system.

(2) Perform the integration of the constitutive model.

(3) Once the integration is done, save the nine components of the rotation matrix (at the end of time increment) in the state variable vector (STATEV).

(4) For the next step, get the rotation matrix (at the beginning of the time increment) from the STATEV vector.

(5) Perform the integration of the constitutive model.

(6)  Update the rotation matrix (i.e. at the end of the time increment) and save it in the STATEV vector.

Repeat (4) to (6) at each time increment.

I am not familiar with the *ORIENTATION option but I don't think it would work since, for crystal plasticity, you need to keep track of the lattice orientation associated with each integration point.

According to my experience, I once implemented VUMAT for transversely isotropic material. I first defined the first direction in VUMAT and then the orientation in the input file (only possible through the input script):

*ORIENTATION *DISTRIBUTION.

Thus, the first direction defined in the subroutine is the orientation of the selected elements in the distribution (in input script) and the second direction (if defined in the subroutine) will make an angle alpha as defined in the subroutine, and so on.

I wish it is helpful!!!

Thanks a lot for your help.

I really do not understand how you defined orientation in subroutine. You can define orientation by *ORIENTATION option but how you interconect it with a variable in the subroutine?

One possible way is to include the three euler angles in your material parameters inside the .inp file. As a result, you will have a set of parameters for each crystal. Then, at the first increment, you calculate the rotation matrix as a function of the euler angles and rotate the single crystal properties in the appropriate coordinate system within your VUMAT.

Thanks Charles, 

Do you mean that I need to take orientation as material property? I thought you cannot update material property! 

What I couldn't infer from Mohammad's answer to my problem is that I really don't know how to use *ORIENTATION option in VUMAT, is there any routine or utility tuat can call these parameters in the subroutine?

Here is how I would do it:

(1) For each integration point, when the VUMAT is called for the first time, read the euler angles (they are stored in the PROPS vector) and calculate the initial rotation matrix (at the beginning of the time increment) that would allow to rotate the single crystal properties in the appropriate coordinate system.

(2) Perform the integration of the constitutive model.

(3) Once the integration is done, save the nine components of the rotation matrix (at the end of time increment) in the state variable vector (STATEV).

(4) For the next step, get the rotation matrix (at the beginning of the time increment) from the STATEV vector.

(5) Perform the integration of the constitutive model.

(6)  Update the rotation matrix (i.e. at the end of the time increment) and save it in the STATEV vector.

Repeat (4) to (6) at each time increment.

I am not familiar with the *ORIENTATION option but I don't think it would work since, for crystal plasticity, you need to keep track of the lattice orientation associated with each integration point.

Thanks Charles; I really appreciate it. I totally understand what is going on. 

What about initializing the STATEV at the first place instead of defining material properties in VUMAT Subroutine? I'm doing polycrystal modelling, so I need to define different material orientation in each material point. I'm going to do this for every material point! 

Technically speaking, using an initialization procedure (such as *INITIAL CONDITIONS) to store the rotation matrix in the STATEV array should work.

However, in my opinion, the most convenient way to account for the initial orientation in a polycrystal consisting of a few hundred crystals is as follows :

- Assign a set of elements to each crystal (elset1=crystal1, elset2=crystal2, ...)

- Define a material (i.e. the PROPOS array) for each crystal such that you have as many material definitions as crystals

Using this procedure, the only difference between 2 different material definitions is the orientation (that is usually specified from 3 euler angles), the other properties should be the same unless your have different metallurgical phases.

Thanks alot for the comments.

HI！ I do not konw the vumat  but in umat there PROPS   The directions in global system and directions in local cubic crystal system of  two nonparallel vectors are needed  to determine the crystal orientation.  so can you send emails  to me the vumat subroutine.   I also want to learning it. thanks a lot  my email  [email protected]

Dear Zhang Jing;

Thanks for your answer to my problem. I'm working on VUMAT code and I'm not finished yet. If you want to code the VUMAT I strongly suggest you to check Prof. Pourboghrat papers on VUMAT.

http://www.sciencedirect.com/science/article/pii/S0749641906001422

The paper is very helpful. I do not need the vumat,    just wnat to the difference between the VUMAT and UMAT  .  If there is no PROPS  about directions  in  the vumat ,you can try to designate the orientation in CAE(inp is better) , I wrote a python program ,which can transform  euler angle into 'inp language'   it's quitre efficient !  [I hope you can understand my Chinglish]

In vumat you just need to update stress tensor. All variables in umat have a substitute here. I think it is a bit simpler here! Because you do not need to calculate elastic_plastic modulus. 

In basic, vumat and umat are different in linearization.

Hi Charles and Ali,

I have written a VUMAT subroutine for single crystal which is based on finite deformation theory . I am facing some problems for Euler angles updates. I am giving Euler angles as initial conditions and from that I am calculating initial rotation matrix.Then integration is done for the constitutive model and  rotation matrix (RN+1=FN+1UN+1-1) is updated for the next time step. Using RN+1,  I am calculating change in Euler angles and for the next step new Euler angles are given as sum of change in Euler angles plus values of Euler angle from previous step. The Euler angles are saved as STATEV vectors.

If I am writing  code like above, I am finding very large change in Euler angles(10 times) even for a simple tension case, which I think is wrong. I really appreciate if you could help me with this problem.

(F=deformation gradient, U=right stretch tensor)

If you have layer dependent material orientation like that of composite laminates, try to use if conditions in your subroutine with the argument "layer".