      subroutine canoq1(ndim,npo,nopoi,coorp,xnorm,np)
c .....................................................
c     IN : 
c     ----
c      ndim : dimension
c      npo  : nombre de points de l'element (3 pour face triangulaire
c                                            4 face quadrangulaire)
c             --- > pour l'instant seul le cas 4 est traite
c             .........................................................
c             pout triangle le calcul de nx,ny,nz est sorti de la boucle
c             et i=1 iq(1,i) = 2 iq(2,i) = 3
c             ..........................................................
c      nopoi(npo)   : numeros des points de l'element (comme ds maillage)
c      coorp(npo,3) : coordonnees 
c      np           : nombre total de points du maillage
c    OUT :
c    -----
c     xnorm(1:3,n)  : somme(n(i)) au point n * aire
c ......................................................
      double precision xnorm(3,np),coorp(npo,ndim),
     *                 x21,y21,x31,y31,z21,z31,delta,nx,ny,nz
      dimension nopoi(npo),iq(2,4)
      data iq /2,4, 3,1, 4,2, 1,3/

      do 1 i = 1,npo

c        n numero global du point
 
         n = nopoi(i)

c        produit vectoriel des aretes adjacentes de i ds le bon sens!       

         X21 = coorp(iq(1,i),1) - coorp(i,1)
         Y21 = coorp(iq(1,i),2) - coorp(i,2)
         X31 = coorp(iq(2,i),1) - coorp(i,1)
         Y31 = coorp(iq(2,i),2) - coorp(i,2)
         Z21 = coorp(iq(1,i),3) - coorp(i,3)
         z31 = coorp(iq(2,i),3) - coorp(i,3)

         nx  =  y21*z31 - z21*y31
         ny  =  z21*x31 - x21*z31
         nz  =  x21*y31 - x31*y21

         xnorm(1,n) = xnorm(1,n)+ nx
         xnorm(2,n) = xnorm(2,n)+ ny
         xnorm(3,n) = xnorm(3,n)+ nz

 1    continue
      end
