c$Id: cubicinpt.f,v 1.2 2002/09/09 23:31:27 sanjay Exp $
      subroutine cubicinpt(ud,nh1,nh3,umat)

c      * * F E A P * * A Finite Element Analysis Program

c....  Copyright (c) 1999: Sanjay Govindjee

c-----[--.----+----.----+----.-----------------------------------------]
c     Purpose: User Constitutive Model Input Routine for a single
c              crystal cubic material

c     Input:
c          umat = 1  
c          umat = 2  

c     Output:
c          ud  [ alpha, beta, gamma, a(3), b(3), c(3) ]
c          nh1     -  Number of gauss point history variables
c          nh3     -  Number of element history variables

c-----[--.----+----.----+----.-----------------------------------------]

      implicit none

      logical pcomp, errck, tinput

      integer nh1,nh3, umat
      character*15 text(2)
      real*8 td(10)

      real*8 aaxis(3),baxis(3),caxis(3),norm

      real*8   ud(*)

      include 'iofile.h'

      save

      text(1) = 'xxxx'
      do while (.not.pcomp(text(1),'ucon',4))
        errck  = tinput(text,1,td,10)

        if (pcomp(text(1),'alpha',5)) then 
            ud(1) = td(1)
        elseif (pcomp(text(1),'beta',4)) then 
            ud(2) = td(1)
        elseif (pcomp(text(1),'gamma',5)) then 
            ud(3) = td(1)
        elseif (pcomp(text(1),'aaxis',5)) then 
            aaxis(1) = td(1)
            aaxis(2) = td(2)
            aaxis(3) = td(3)
        elseif (pcomp(text(1),'baxis',5)) then 
            baxis(1) = td(1)
            baxis(2) = td(2)
            baxis(3) = td(3)
        end if

      end do


      norm = dsqrt(aaxis(1)**2.d0+aaxis(2)**2.d0+aaxis(3)**2.d0)
      ud(4) = aaxis(1)/norm
      ud(5) = aaxis(2)/norm
      ud(6) = aaxis(3)/norm
      norm = dsqrt(baxis(1)**2.d0+baxis(2)**2.d0+baxis(3)**2.d0)
      ud(7) = baxis(1)/norm
      ud(8) = baxis(2)/norm
      ud(9) = baxis(3)/norm
      caxis(1) = aaxis(2)*baxis(3)-aaxis(3)*baxis(2)
      caxis(2) = -aaxis(1)*baxis(3)+aaxis(3)*baxis(1)
      caxis(3) = aaxis(1)*baxis(2)-aaxis(2)*baxis(1)
      norm = dsqrt(caxis(1)**2.d0+caxis(2)**2.d0+caxis(3)**2.d0)
      ud(10) = caxis(1)/norm
      ud(11) = caxis(2)/norm
      ud(12) = caxis(3)/norm

     
      nh1 = 0
      nh3 = 0

      if(ior.lt.0) then
        write(*,2001) 
        write(*,2003) ud(1),ud(2),ud(3)
        write(*,2004) ud(4),ud(5),ud(6)
        write(*,2005) ud(7),ud(8),ud(9)
        write(*,2006) ud(10),ud(11),ud(12)
      else
        write(iow,2001) 
        write(iow,2003) ud(1),ud(2),ud(3)
        write(iow,2004) ud(4),ud(5),ud(6)
        write(iow,2005) ud(7),ud(8),ud(9)
        write(iow,2006) ud(10),ud(11),ud(12)
      endif

2001  format(/5x,'Single Crystal Cubic Material'/)
2003  format(/5x,'Alpha  ',e12.5/
     &        5x,'Beta   ',e12.5/
     &        5x,'Gamma  ',e12.5)
2004  format(/5x,'A-Axis ',3e14.3/)
2005  format(5x,'B-Axis ',3e14.3/)
2006  format(5x,'C-Axis ',3e14.3/)

      end