src/tool/Utilitaire/uttoqu.F

   1       subroutine uttoqu (s1, s2, s3, s4, coonoe, torsio)
   2 c ______________________________________________________________________
   3 c
   4 c                             H O M A R D
   5 c
   6 c Outil de Maillage Adaptatif par Raffinement et Deraffinement d'EDF R&D
   7 c
   8 c Version originale enregistree le 18 juin 1996 sous le numero 96036
   9 c aupres des huissiers de justice Simart et Lavoir a Clamart
  10 c Version 11.2 enregistree le 13 fevrier 2015 sous le numero 2015/014
  11 c aupres des huissiers de justice
  12 c Lavoir, Silinski & Cherqui-Abrahmi a Clamart
  13 c
  14 c    HOMARD est une marque deposee d'Electricite de France
  15 c
  16 c Copyright EDF 1996
  17 c Copyright EDF 1998
  18 c Copyright EDF 2002
  19 c Copyright EDF 2020
  20 c ______________________________________________________________________
  21 c
  22 c     UTilitaire : TOrsion d'un QUadrangle
  23 c     --           --           --
  24 c ______________________________________________________________________
  25 c .        .     .        .                                            .
  26 c .  nom   . e/s . taille .           description                      .
  27 c .____________________________________________________________________.
  28 c .s1/2/3/4. e   .  1     . les noeuds de la face                      .
  29 c . coonoe . e   . nbnoto . coordonnees des noeuds                     .
  30 c .        .     . * sdim .                                            .
  31 c . torsio .  s  .  1     . torsion de la face entre 0, plane, et 2    .
  32 c .____________________________________________________________________.
  33 c
  34 c====
  35 c 0. declarations et dimensionnement
  36 c====
  37 c
  38 c 0.1. ==> generalites
  39 c
  40       implicit none
  41       save
  42 c
  43 c 0.2. ==> communs
  44 c
  45 #include "nombno.h"
  46 c
  47 c 0.3. ==> arguments
  48 c
  49       integer s1, s2, s3, s4
  50 c
  51       double precision coonoe(nbnoto,3)
  52       double precision torsio
  53 c
  54 c 0.4. ==> variables locales
  55 c
  56       integer iaux
  57 c
  58       double precision n1(3), n2(3)
  59       double precision vs1s2(3), vs1s3(3)
  60       double precision vs1s4(3), vs2s4(3), vs3s4(3)
  61       double precision daux
  62 c
  63 #include "impr03.h"
  64 c
  65 c====
  66 c 1. Traitement par la diagonale s1/s3
  67 c====
  68 c
  69 c          s1                       s2
  70 c           -------------------------
  71 c           .   .                   .
  72 c           .       .               .
  73 c           .           .           .
  74 c           .               .       .
  75 c           .                   .   .
  76 c           -------------------------
  77 c          s4                       s3
  78 c   Attention a prendre la meme orientation pour les deux normales !
  79 c
  80 c 1.1. ==> Vecteur normal au triangle (s1,s2,s3)
  81 c
  82 cgn      write(*,90002) 'triangle des sommets', s1, s2, s3
  83       vs1s2(1) = coonoe(s2,1) - coonoe(s1,1)
  84       vs1s2(2) = coonoe(s2,2) - coonoe(s1,2)
  85       vs1s2(3) = coonoe(s2,3) - coonoe(s1,3)
  86 c
  87       vs1s3(1) = coonoe(s3,1) - coonoe(s1,1)
  88       vs1s3(2) = coonoe(s3,2) - coonoe(s1,2)
  89       vs1s3(3) = coonoe(s3,3) - coonoe(s1,3)
  90 c
  91       n1(1) = vs1s3(2)*vs1s2(3) - vs1s3(3)*vs1s2(2)
  92       n1(2) = vs1s3(3)*vs1s2(1) - vs1s3(1)*vs1s2(3)
  93       n1(3) = vs1s3(1)*vs1s2(2) - vs1s3(2)*vs1s2(1)
  94 c
  95       daux = sqrt(n1(1)*n1(1)+n1(2)*n1(2)+n1(3)*n1(3))
  96       do 11 , iaux = 1 , 3
  97         n1(iaux) = n1(iaux)/daux
  98    11 continue
  99 cgn      write(*,92010) 11, n1
 100 c
 101 c 1.2. ==> Vecteur normal au triangle (s1,s4,s3)
 102 c
 103 cgn      write(*,90002) 'triangle des sommets', s1, s3, s4
 104       vs1s4(1) = coonoe(s4,1) - coonoe(s1,1)
 105       vs1s4(2) = coonoe(s4,2) - coonoe(s1,2)
 106       vs1s4(3) = coonoe(s4,3) - coonoe(s1,3)
 107 c
 108       n2(1) = vs1s3(2)*vs1s4(3) - vs1s3(3)*vs1s4(2)
 109       n2(2) = vs1s3(3)*vs1s4(1) - vs1s3(1)*vs1s4(3)
 110       n2(3) = vs1s3(1)*vs1s4(2) - vs1s3(2)*vs1s4(1)
 111 c
 112       daux = sqrt(n2(1)*n2(1)+n2(2)*n2(2)+n2(3)*n2(3))
 113       do 12 , iaux = 1 , 3
 114         n2(iaux) = -n2(iaux)/daux
 115    12 continue
 116 cgn      write(*,92010) 12, n2
 117 c
 118 c 1.3. ==> Produit scalaire des vecteurs normaux
 119 c
 120       torsio = n1(1)*n2(1) + n1(2)*n2(2) + n1(3)*n2(3)
 121 cgn      write(*,92010) n1(1)*n2(1) + n1(2)*n2(2) + n1(3)*n2(3)
 122 c
 123 c====
 124 c 2. Traitement par la diagonale s2/s4
 125 c====
 126 c
 127 c          s1                       s2
 128 c           -------------------------
 129 c           .                   .   .
 130 c           .               .       .
 131 c           .           .           .
 132 c           .       .               .
 133 c           .   .                   .
 134 c           -------------------------
 135 c          s4                       s3
 136 c   Attention a prendre la meme orientation pour les deux normales !
 137 c
 138 c 2.1. ==> Vecteur normal au triangle (s4,s1,s2)
 139 c
 140 cgn      write(*,90002) 'triangle des sommets', s4,s1,s2
 141       vs2s4(1) = coonoe(s4,1) - coonoe(s2,1)
 142       vs2s4(2) = coonoe(s4,2) - coonoe(s2,2)
 143       vs2s4(3) = coonoe(s4,3) - coonoe(s2,3)
 144 c
 145       n1(1) = vs2s4(2)*vs1s4(3) - vs2s4(3)*vs1s4(2)
 146       n1(2) = vs2s4(3)*vs1s4(1) - vs2s4(1)*vs1s4(3)
 147       n1(3) = vs2s4(1)*vs1s4(2) - vs2s4(2)*vs1s4(1)
 148 c
 149       daux = sqrt(n1(1)*n1(1)+n1(2)*n1(2)+n1(3)*n1(3))
 150       do 21 , iaux = 1 , 3
 151         n1(iaux) = n1(iaux)/daux
 152    21 continue
 153 cgn      write(*,92010) 21, n1
 154 c
 155 c 2.2. ==> Vecteur normal au triangle (s1,s4,s3)
 156 c
 157 cgn      write(*,90002) 'triangle des sommets', s4, s2, s3
 158       vs3s4(1) = coonoe(s4,1) - coonoe(s3,1)
 159       vs3s4(2) = coonoe(s4,2) - coonoe(s3,2)
 160       vs3s4(3) = coonoe(s4,3) - coonoe(s3,3)
 161 c
 162       n2(1) = vs3s4(2)*vs2s4(3) - vs3s4(3)*vs2s4(2)
 163       n2(2) = vs3s4(3)*vs2s4(1) - vs3s4(1)*vs2s4(3)
 164       n2(3) = vs3s4(1)*vs2s4(2) - vs3s4(2)*vs2s4(1)
 165 c
 166       daux = sqrt(n2(1)*n2(1)+n2(2)*n2(2)+n2(3)*n2(3))
 167       do 22 , iaux = 1 , 3
 168         n2(iaux) = n2(iaux)/daux
 169    22 continue
 170 cgn      write(*,92010) 22, n2
 171 c
 172 c 2.3. ==> Produit scalaire des vecteurs normaux
 173 c
 174       torsio = min ( torsio, n1(1)*n2(1) + n1(2)*n2(2) + n1(3)*n2(3) )
 175 cgn      write(*,92010) n1(1)*n2(1) + n1(2)*n2(2) + n1(3)*n2(3)
 176 c
 177 c====
 178 c 3. Torsion : Le produit scalaire des vecteurs normaux vaut 1 s'ils
 179 c    sont colineaires. La torsion vaudra alors 0. A l'extreme, pour
 180 c    deux triangles opposes, le produit scalaire vaut -1 et la
 181 c    torsion vaudra 2.
 182 c    Remarque : on prend la valeur absolue pour eviter les affichages
 183 c               en -0.00000
 184 c====
 185 c
 186       torsio = abs(1.d0 - torsio)
 187 c
 188       end