sortho31.F File Reference

#include "implicit_f.inc"

Go to the source code of this file.

Functions/Subroutines
subroutine	sortho31 (x1, x2, x3, x4, x5, x6, x7, x8, y1, y2, y3, y4, y5, y6, y7, y8, z1, z2, z3, z4, z5, z6, z7, z8, e1x, e2x, e3x, e1y, e2y, e3y, e1z, e2z, e3z)

Function/Subroutine Documentation

sortho31()

subroutine sortho31	(		x1,
			x2,
			x3,
			x4,
			x5,
			x6,
			x7,
			x8,
			y1,
			y2,
			y3,
			y4,
			y5,
			y6,
			y7,
			y8,
			z1,
			z2,
			z3,
			z4,
			z5,
			z6,
			z7,
			z8,
			e1x,
			e2x,
			e3x,
			e1y,
			e2y,
			e3y,
			e1z,
			e2z,
			e3z )

Definition at line 29 of file sortho31.F.

C-----------------------------------------------                        
C   I m p l i c i t   T y p e s                                         
C-----------------------------------------------                        
#include      "implicit_f.inc"                                          
C-----------------------------------------------                        
C   D u m m y   A r g u m e n t s                                       
C-----------------------------------------------                        
      my_real                                                           
     .   x1, x2, x3, x4, x5, x6, x7, x8,                                
     .   y1, y2, y3, y4, y5, y6, y7, y8,                                
     .   z1, z2, z3, z4, z5, z6, z7, z8,                                
     .   e1x, e1y, e1z,                                                 
     .   e2x, e2y, e2z,                                                 
     .   e3x, e3y, e3z                                                  
C-----------------------------------------------                        
C   L o c a l   V a r i a b l e s                                       
C-----------------------------------------------                        
      INTEGER N,NITER                                                   
      my_real                                                           
     .   x17 , x28 , x35 , x46,                                         
     .   y17 , y28 , y35 , y46,                                         
     .   z17 , z28 , z35 , z46,                                         
     .   a17 , a28 ,                                                    
     .   b17 , b28 ,                                                    
     .   c17 , c28 ,                                                    
     .   rx , ry , rz ,                                                 
     .   sx , sy , sz ,                                                 
     .   tx , ty , tz                                                   
      my_real                                                           
     .   aa,bb                                                          
      DATA niter/3/                                                     
C-----------------------------------------------                        
      x17=x7-x1                                                         
      x28=x8-x2                                                         
      x35=x5-x3                                                         
      x46=x6-x4                                                         
      y17=y7-y1                                                         
      y28=y8-y2                                                         
      y35=y5-y3                                                         
      y46=y6-y4                                                         
      z17=z7-z1                                                         
      z28=z8-z2                                                         
      z35=z5-z3                                                         
      z46=z6-z4                                                         
      rx=x17+x28-x35-x46                                                
      ry=y17+y28-y35-y46                                                
      rz=z17+z28-z35-z46                                                
      a17=x17+x46                                                       
      a28=x28+x35                                                       
      b17=y17+y46                                                       
      b28=y28+y35                                                       
      c17=z17+z46                                                       
      c28=z28+z35                                                       
      sx=a17+a28                                                        
      sy=b17+b28                                                        
      sz=c17+c28                                                        
      tx=a17-a28                                                        
      ty=b17-b28                                                        
      tz=c17-c28                                                        
c~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~                       
c     norme r s t                                                       
c~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~                       
      aa = sqrt(rx*rx + ry*ry + rz*rz)                                  
      if ( aa/=zero) aa = one / aa                                     
      rx = rx * aa                                                      
      ry = ry * aa                                                      
      rz = rz * aa                                                      
      aa = sqrt(sx*sx + sy*sy + sz*sz)                                  
      if ( aa/=zero) aa = one / aa                                     
      sx = sx * aa                                                      
      sy = sy * aa                                                      
      sz = sz * aa                                                      
      aa = sqrt(tx*tx + ty*ty + tz*tz)                                  
      if ( aa/=zero) aa = one / aa                                     
      tx = tx * aa                                                      
      ty = ty * aa                                                      
      tz = tz * aa                                                      
c~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~                       
c     iterations                                                        
c~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~                       
      n=0                                                               
111   CONTINUE                                                          
      n=n+1                                                             
c~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~                       
       e1x = sy * tz - sz * ty + rx                                     
       e1y = sz * tx - sx * tz + ry                                     
       e1z = sx * ty - sy * tx + rz                                     
c                                                                       
       e2x = ty * rz - tz * ry + sx                                     
       e2y = tz * rx - tx * rz + sy                                     
       e2z = tx * ry - ty * rx + sz                                     
c                                                                       
       e3x = ry * sz - rz * sy + tx                                     
       e3y = rz * sx - rx * sz + ty                                     
       e3z = rx * sy - ry * sx + tz                                     
c                                                                       
       bb = sqrt(e1x*e1x + e1y*e1y + e1z*e1z)                           
       if ( bb/=zero) bb = one / bb                                    
       rx = e1x * bb                                                    
       ry = e1y * bb                                                    
       rz = e1z * bb                                                    
c                                                                       
       bb = sqrt(e2x*e2x + e2y*e2y + e2z*e2z)                           
       if ( bb/=zero) bb = one / bb                                    
       sx = e2x * bb                                                    
       sy = e2y * bb                                                    
       sz = e2z * bb                                                    
c                                                                       
       bb = sqrt(e3x*e3x + e3y*e3y + e3z*e3z)                           
       if ( bb/=zero) bb = one / bb                                    
       tx = e3x * bb                                                    
       ty = e3y * bb                                                    
       tz = e3z * bb                                                    
c                                                                       
      IF (n<niter) GOTO 111                                          
c     norme et orthogonalisation                                        
c~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~                       
      e1x = rx                                                          
      e1y = ry                                                          
      e1z = rz                                                          
c                                                                       
      e3x = e1y * sz - e1z * sy                                         
      e3y = e1z * sx - e1x * sz                                         
      e3z = e1x * sy - e1y * sx                                         
c                                                                       
      aa = sqrt(e3x*e3x + e3y*e3y + e3z*e3z)                            
      if ( aa/=zero) aa = one / aa                                     
      e3x = e3x * aa                                                    
      e3y = e3y * aa                                                    
      e3z = e3z * aa                                                    
c                                                                       
      e2x = e3y * e1z - e3z * e1y                                       
      e2y = e3z * e1x - e3x * e1z                                       
      e2z = e3x * e1y - e3y * e1x                                       
c~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~                       
      RETURN