c     mkFem: start of new addmel.f 
c     
      subroutine AMEDSP(ivar, elnum, nodpek, nnode, mnnode,
     $             coord, ntunk, u, uo, atw, bindx, idxatw)
      include 'include/size.h'
      include 'slprms.h' 
      include 'eqprms.h' 
  
      integer ivar, elnum, nnode, mnnode, ntunk
      integer bindx(0:*) 
      double precision atw(0:*) 
      double precision u (0:ntunk-1,nvar)
      double precision uo(0:ntunk-1,nvar)
      double precision coord(0:mnnode-1,spcdim)
      integer nodpek(elnum,nobf,noofmx)
      integer idxatw(elnum, nobf, nobf)
  
c     local integer variables   
      integer ind    
  
      do ind=0,bindx(nnode)-1  
        atw(ind) = 0.d0  
      enddo  
      if (ivar.eq.1) then
          call addm1(elnum,nodpek,nnode,mnnode,coord,ntunk,
     $              u,uo,atw,bindx,idxatw)
      endif 
      return 
      end 
c########################## 
      subroutine adidsp(ivar,dinv)  
  
      print *, 'adidsp: error: adidsp should not be called '  
        
      return 
      end 
 
      subroutine UZADSP(ivar, elnum, nodpek, nnode, mnnode,
     $             coord, ntunk, u, uo, atw, bindx, idxatw)
      include 'include/size.h'
      include 'slprms.h' 
      include 'eqprms.h' 
  
      integer ivar, elnum, nnode, mnnode, ntunk
      integer bindx(0:*) 
      double precision atw(0:*) 
      double precision u (0:ntunk-1,nvar)
      double precision uo(0:ntunk-1,nvar)
      double precision coord(0:mnnode-1,spcdim)
      integer nodpek(elnum,nobf,noofmx)
      integer idxatw(elnum, nobf, nobf)
  
c     local integer variables   
      integer ind    
  
      do ind=0,bindx(nnode)-1  
        atw(ind) = 0.d0  
      enddo  
      return 
      end 
c########################## 
 
C ####################################### 
      subroutine addm1(elnum,nodpek,nnode,mnnode,
     $     coord,ntunk,u,uo,atw,bindx,idxatw)
      implicit double precision (t)  
      implicit double precision (s)  
      include 'include/size.h'
      include 'slprms.h' 
      include 'eqprms.h' 
      include 'grid.h' 
      
      integer elnum, nnode, ntunk, mnnode
      integer bindx(0:*)
      double precision atw(0:*)
      double precision u (0:mnnode-1, nvar)
      double precision uo(0:mnnode-1, nvar)
      double precision coord(0:mnnode-1,spcdim)
      integer          nodpek(elnum,nobf,noofmx)
      integer          idxatw(elnum,nobf,nobf)
 
c     local integer variables   
      double precision Atmp(nobf,nobf)  
      integer el,li,lj    
      integer t1,t3,t6
      ivar=1
      do el =1,elnum
      t1 = nodpek(el,3,1)
      t2 = coord(t1,1)
      t3 = nodpek(el,2,1)
      t4 = coord(t3,2)
      t6 = nodpek(el,1,1)
      t7 = coord(t6,2)
      t9 = coord(t6,1)
      t11 = coord(t3,1)
      t12 = coord(t1,2)
      t16 = -t2*t4+t2*t7+t9*t4+t11*t12-t11*t7-t9*t12
      t18 = dabs(t16)/2
      t19 = -t4+t12
      t20 = t18*t19
      t21 = t16**2
      t22 = 1/t21
      t24 = -t22*t19
      t26 = -t11+t2
      t27 = t18*t26
      t29 = -t22*t26
      Atmp(1,1) = -t20*t24-t27*t29
      t31 = -t12+t7
      t32 = t18*t31
      t34 = -t2+t9
      t35 = t18*t34
      Atmp(1,2) = -t32*t24-t35*t29
      t37 = t4-t7
      t38 = t18*t37
      t40 = t11-t9
      t41 = t18*t40
      Atmp(1,3) = -t38*t24-t41*t29
      t43 = t31*t22
      t45 = t34*t22
      Atmp(2,1) = t20*t43+t27*t45
      t47 = t31**2
      t50 = t34**2
      Atmp(2,2) = t18*t47*t22+t18*t50*t22
      Atmp(2,3) = t38*t43+t41*t45
      t56 = -t37*t22
      t59 = -t40*t22
      Atmp(3,1) = -t20*t56-t27*t59
      Atmp(3,2) = -t32*t56-t35*t59
      Atmp(3,3) = -t38*t56-t41*t59
        do li=1,3 
c       If the nodes are in this processor 
           if (nodpek(el,li,1).lt.nnode) then 
             do lj=1,3 
                 atw(idxatw(el,li,lj))=atw(idxatw(el,li,lj))+Atmp(li,lj)
              enddo 
           endif 
        enddo 
      enddo 
      return  
      end  
