/* c$Header: /usr/home/sanjay/Feap/Elem/Ec/RCS/stcn13.c,v 1.2 2000/01/28 20:42:11 sanjay Exp $ */ 

#include <stdlib.h>

/* zero index version of st */
#define st(a,b) (*(stf+a+(b)*numnp))


void stcn13(int *ix, double shps[4][4][3], double *xsj,
            double sigl[4][6],
            double *dt, double *stf,
            int lint,int nel, int numnp) {

/*
c-------------------------------------------------------------------
c
c      * * F E A P * * A Finite Element Analysis Program
c      Copyright (c) 1984-1994: Robert L. Taylor
c
c      C Language Example Element for FEAP
c      Copyright (c) 2000     : Sanjay Govindjee
c
c      Stress projection routine
c
c      Inputs:
c         ix               - connectivity information
c         shps(lint,nel,3) - shape fncs and derivatives
c         xsj(lint)        - Gauss point jacobians
c         sigl(lint,6)     - Gauss point stresses
c         lint             - number of Gauss points
c         nel              - number of element nodes
c         numnp            - number of nodes in problem
c
c      Output:
c         dt(numnp)        - nodal weights
c         st(numnp,6)      - nodal stresses (accessed throught cover
c                            function and pointer stf)
c
c-------------------------------------------------------------------
*/

      int i, l, ll;

      double  xj;

      /* compute stress projections to nodes */
      for(l = 0 ; l < lint ; l++ ) {

      /* compute projections: int ( sig * shp(i) * darea ) */

        for(i = 0 ; i < nel ; i++ ) {
          ll = abs(ix[i]);
          if(ll != 0) {
            ll -= 1;   /* account for zero based indexing in C */
            xj       = xsj[l]*shps[l][i][2];
            dt[ll]   +=  xj;
            st(ll,0) +=  sigl[l][0]*xj;
            st(ll,1) +=  sigl[l][1]*xj;
            st(ll,2) +=  sigl[l][2]*xj;
            st(ll,3) +=  sigl[l][3]*xj;
            st(ll,4) +=  sigl[l][4]*xj;
            st(ll,5) +=  sigl[l][5]*xj;
          }
        }
      }


}