User element with 5 dof gives wrong output

[pradeep]

Hi,
I will take you through a sequence of things I did for writing a user element with 5dof.
1. I wrote a user element with 3 dof/node(solid element with ux,uy,uz being the dof) and this worked properly.
2. Then I wrote a user element with 1 dof/node which is temperature dof and this worked as expected as well.
3. Next I wrote a user element with 1dof/ node for normalized concentration and this owrked as well.
4. Now i wanted to combine all the 5 dofs togather .So I coupled solid and thermal elements first which worked without any hiccup.Next when i tried to add diffusion element even the thermal dof stopped working.To debug the problem I sent Null values of residual and mass matrix corresponding to diffusion before the assembly of elemental stiffness matrix.Still the temperature dof doesnot work.
First 5 rows of my stiffness matrix would look like this after assemble:
UX1 UY1 UZ1 0 0 UX2 UY2 UZ2 0 0................UX8 UY8 UZ8 0 0
UX1 UY1 UZ1 0 0 UX2 UY2 UZ2 0 0................UX8 UY8 UZ8 0 0
UX1 UY1 UZ1 0 0 UX2 UY2 UZ2 0 0................UX8 UY8 UZ8 0 0
0       0        0     T1  0         0       0   T2 0................0       0         0 T8 0
0       0        0     0    D1         0       0   0 D2................0       0        0 0  D8
First row of my residual vector after assembly
RUX1 RUY1 RUZ1 RT1 RD1.................................RUX8 RUY8 RUZ8 RT8 RD8
I am writing it for a 8 node brick element and for a  weakly coupled system.
I am attaching my input file as well.
Thanks in advance.
     

[pradeep]

Another addition to the problem is that even if I give the following boundary conditions i should get zero values for temperature degree of freedom as it is a weakly coupled system.But I get values slowly going unconverged.I am attachina a snapshop of Temperature contour plot .
EBOU
  1 0 1 0 0 0 0
  1 10 1 1 1 1 0
   
EDIS
  1 0 1 0 0 0 0

[Prof. R.L. Taylor]

There does not seem to be any reason why 5 dof does not work.  FEAPpv can handle shells and beams which have this many.

First: What data does the output file give?  Is everything there you expect?

Second: what happens if you run your solid element with 5 dof specified on the control record.

Assuming all this works.  You need to start output for element results.  Do a 1- elmeent problem until you find what is going wrong.  Note that utilities "mprint" and "iprint" can help output of real and integer results.

Are you handling the 'ndf x ndf' blocks correctly in defining the tangent and residual?

[pradeep]

Thank you professor Taylor.
I restructured the code ,initializing the variables properly and the code is working now.

As I was using double precision matrices , identity matrix minus deformation gradient (I-F) was giving me 1x10^-16 values.So i wrote a subroutine to reduce the precision.

May be i was going wrong with multiplication of matrices somewhere and restructuring it again helped.

Thank You!

[Prof. R.L. Taylor]

Glad you found the problem.  In the latest release for constitution we compute H = F - I directly from the displacement state for as many problems as possible, thus avoiding the roundoff problem.  In our computation of kinematics we compute F(3,3,4): F(:,:,1) = F_n+1; F(:,:,2) = F_n; F(:,:,3) = H_n+1; F(:,:,4) = H_n -- same for detF.  If you have latest ver8.4 you can look in the elements to find the computations. 

[pradeep]

Thank You Professor,
I will have a look.