daily pastebin goal
1%
SHARE
TWEET

Untitled

nirn May 16th, 2018 124 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. #include <stdio.h>
  2. #include <stdlib.h>
  3.  
  4. struct vector
  5. {
  6.     float data[4];
  7. };
  8.  
  9. float sumVec(struct vector v){
  10.     return v.data[0]+v.data[1]+v.data[2]+v.data[3];
  11. }
  12.  
  13. float det(struct vector matrix[4])
  14. {  
  15.     // matrix - tablica czterech wektorów w komentarzach kolejno A, B, C i D
  16.     struct vector m;
  17.     float result;
  18.  
  19.     /* Obliczenia do wykonania:
  20.     det =   (a0*b1-a1*b0)*(c2*d3-c3*d2) +
  21.         (a2*b0-a0*b2)*(c1*d3-c3*d1) +
  22.         (a0*b3-a3*b0)*(c1*d2-c2*d1) +
  23.         (a1*b2-a2*b1)*(c0*d3-c3*d0) +
  24.  
  25.         (a2*b3-a3*b2)*(c0*d1-c1*d0) +
  26.         (a3*b1-a1*b3)*(c0*d2-c2*d0)
  27.  
  28.     Jak widać w druga suma jest zbyt mała by uzupełnić nią 4-elementowy wektor
  29.     Stosujemy zatem następujące przekształcenia:
  30.  
  31.     det =   (a0*b1-a1*b0)*(c2*d3-c3*d2) +
  32.         (a2*b0-a0*b2)*(c1*d3-c3*d1) +
  33.         (a0*b3-a3*b0)*(c1*d2-c2*d1) +
  34.         (a1*b2-a2*b1)*(c0*d3-c3*d0) +
  35.  
  36.         (a2*b3-a3*b2)*(c0*d1-c1*d0) +
  37.         (a3*b1-a1*b3)*(c0*d2-c2*d0) +
  38.               0      *      0       +
  39.               0      *      0
  40.  
  41.     det =   (a0*b1-a1*b0)*(c2*d3-c3*d2) +
  42.         (a2*b0-a0*b2)*(c1*d3-c3*d1) +
  43.         (a0*b3-a3*b0)*(c1*d2-c2*d1) +
  44.         (a1*b2-a2*b1)*(c0*d3-c3*d0) +
  45.  
  46.         (a2*b3-a3*b2)*(c0*d1-c1*d0) +
  47.         (a3*b1-a1*b3)*(c0*d2-c2*d0) +
  48.         (a0*b0-a0*b0)*(c0*d0-c0*d0) +
  49.         (a0*b0-a0*b0)*(c0*d0-c0*d0)
  50.  
  51.     Teraz kolumny indeksów elementów wektora mogą stać się argumentami kolejnych wywołań shufps.
  52.     Każdy wektor przestawimy zgodnie z indeksami a następnie będziemy wykonywać operacje
  53.     na wektorach czteroelementowych.
  54.     */
  55.    
  56.     asm (                                              
  57.         /* Obliczenie 1 składowego wyniku
  58.         (a0*b1-a1*b0)
  59.         (a2*b0-a0*b2)
  60.         (a0*b3-a3*b0)
  61.         (a1*b2-a2*b1)
  62.          A' B' A" B"
  63.         */
  64.         "movups (%0), %%xmm0     \n\t"
  65.         "movups (%1), %%xmm1 \n\t"  
  66.        
  67.         "shufps $0x48 , %%xmm0, %%xmm0 \n\t"            //A' (1, 0, 2, 0)
  68.         "shufps $0xb1 , %%xmm1, %%xmm1 \n\t"            //B' (2, 3, 0, 1)
  69.         "mulps %%xmm0,%%xmm1 \n\t"              //A'*B'
  70.        
  71.         // ponowne ustawienie zawartości rejestrów
  72.         // Obliczenie iloczynów z prawej strony
  73.         "movups (%0), %%xmm0 \n\t"             
  74.         "movups (%1), %%xmm2 \n\t"             
  75.        
  76.         "shufps $0xb1,%%xmm0,%%xmm0 \n\t"           //A" (2, 3, 0, 1)
  77.         "shufps $0x48,%%xmm2,%%xmm2 \n\t"           //B" (1, 0, 2, 0)
  78.         "mulps  %%xmm0,%%xmm2 \n\t"                 //A"*B"
  79.        
  80.         // Wyliczenie różnicy
  81.         "subps %%xmm1, %%xmm2 \n\t"                 //A'*B'-A"*B" po wektoryzacji
  82.        
  83.         /*
  84.         (c2*d3-c3*d2)
  85.         (c1*d3-c3*d1)
  86.         (c1*d2-c2*d1)
  87.         (c0*d3-c3*d0)
  88.          C' D' C" D"
  89.         */ 
  90.         "movups (%2), %%xmm0     \n\t"              //C,
  91.         "movups (%3), %%xmm1 \n\t"                  //D
  92.        
  93.         "shufps $0x16 , %%xmm0, %%xmm0 \n\t"            //C' (0, 1, 1, 2)
  94.         "shufps $0xef , %%xmm1, %%xmm1 \n\t"            //D' (3, 2, 3, 3)
  95.         "mulps %%xmm0,%%xmm1 \n\t"              // C'*D'
  96.        
  97.         "movups (%2), %%xmm0 \n\t"              //C
  98.         "movups (%3), %%xmm3 \n\t"                  //D
  99.        
  100.         "shufps $0xef,%%xmm0 ,%%xmm0 \n\t"              //C" (3, 2, 3, 3)
  101.         "shufps $0x16, %%xmm3, %%xmm3 \n\t"             //D" (0, 1, 1, 2)
  102.         "mulps  %%xmm0,%%xmm3 \n\t"                 //C"*D"
  103.        
  104.         "subps %%xmm1, %%xmm3 \n\t"                 //C'*D'-C"*D"
  105.  
  106.         "mulps %%xmm2,%%xmm3 \n\t"              //(A'*B'-A"*B")*(C'*D'-C"*D")  - wynik skladowy 1
  107.    
  108.    
  109.         /* Obliczenie 2 składowego wyniku
  110.         (a2*b3-a3*b2)
  111.         (a3*b1-a1*b3)
  112.         (a0*b0-a0*b0) <- 0
  113.         (a0*b0-a0*b0) <- 0
  114.          A' B' A" B"
  115.         */
  116.         "movups (%0), %%xmm0     \n\t"              //A
  117.         "shufps $0xe , %%xmm0, %%xmm0  \n\t"            //A' (0, 0, 3, 2)
  118.        
  119.         "movups (%1), %%xmm2 \n\t"                  //B
  120.         "shufps $0x7 , %%xmm2, %%xmm2 \n\t"             //B' (0, 0, 1, 3)
  121.        
  122.         "mulps %%xmm0, %%xmm2 \n\t"                 //A'*B'
  123.        
  124.         "movups (%0), %%xmm0     \n\t"              //A
  125.         "shufps $0x7 , %%xmm0, %%xmm0 \n\t"             //A" (0, 0, 1, 3)
  126.        
  127.         "movups (%1), %%xmm4     \n\t"              //B
  128.         "shufps $0xe , %%xmm4, %%xmm4 \n\t"             //B" (0, 0, 3, 2)
  129.        
  130.         "mulps %%xmm0,%%xmm4 \n\t"              //A"*B"
  131.        
  132.         "subps %%xmm4, %%xmm2 \n\t"                 //A'*B'-A"*B"
  133.        
  134.         /*
  135.         (c0*d1-c1*d0)
  136.         (c0*d2-c2*d0)
  137.         (c0*d0-c0*d0) <- 0
  138.         (c0*d0-c0*d0) <- 0
  139.          C' D' C" D"
  140.         */
  141.         "movups (%2), %%xmm0 \n\t"              //C
  142.         "shufps $0x0 , %%xmm0, %%xmm0 \n\t"         //C' (0, 0, 0, 0)
  143.        
  144.         "movups (%3), %%xmm4 \n\t"                  //D
  145.         "shufps $0x9 , %%xmm4, %%xmm4 \n\t"             //D' (0, 0, 2, 1)
  146.        
  147.         "mulps %%xmm0, %%xmm4 \n\t"                 //C * D
  148.        
  149.         "movups (%2), %%xmm0 \n\t"              //C
  150.         "shufps $0x9 , %%xmm0, %%xmm0 \n\t"         //C" (0, 0, 2, 1)
  151.        
  152.         "movups (%3), %%xmm5 \n\t"                  //D
  153.         "shufps $0x0 , %%xmm5, %%xmm5 \n\t"         //D" (0, 0, 0, 0)
  154.        
  155.         "mulps %%xmm5,%%xmm0 \n\t"              //C" * D"
  156.        
  157.         "subps %%xmm4, %%xmm0 \n\t"                 //C'*D'- C"*D"
  158.        
  159.         "mulps %%xmm0,%%xmm2 \n\t"              //skladowa wyniku 2
  160.        
  161.         "subps %%xmm2,%%xmm3 \n\t"              //roznica skladowych = wynik
  162.         "movups %%xmm3,(%4) \n\t"               //wektor rozwiazan
  163.         :
  164.         : "r" (matrix[0].data),"r" (matrix[1].data) , "r" (matrix[2].data), "r"(matrix[3].data), "r" (m.data)
  165.     );
  166.    
  167.     int i;
  168.     printf("Sumy częściowe: \n");
  169.     for(i = 0; i<4; i++)
  170.         printf("%f \n",m.data[i]);
  171.  
  172.     return sumVec(m);
  173. };
  174.  
  175.  
  176. int main()
  177. {
  178.  
  179.     float deter;
  180.  
  181.     struct vector matrix[4];
  182.  
  183.     // tablica użyta jedynie do ułatwienia zapisu macierzowego
  184.     // matrix można uzupełniać przez matrix[i].data[j] bez inicjalizacji tablicy
  185.     float m[4][4] = {{1, 2, 3, 4},
  186.              {5, 6, 7, 8},
  187.              {9, 1, 2, 3},
  188.              {4, 5, 6, 8} };   
  189.  
  190.     int i, j;
  191.  
  192.     for(i = 0; i < 4; i++){
  193.         for(j = 0; j < 4;  j++)
  194.             matrix[i].data[j] = m[i][j];
  195.     }
  196.    
  197.  
  198.     deter=det(matrix);
  199.     printf("Wyznacznik: \n");
  200.     printf("%f \n",deter); 
  201.    
  202.     return 0;
  203. }
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy. OK, I Understand
 
Top